I want to take advantage of this small space to thank several people whose help has been essential for completing my thesis. First of all, I want to thank my supervisors Pol Bernard Gossiaux and Joerg Aichelin, for their advice, contributions, corrections, comments, and contradictory criticisms that have helped me evolve as a scientist and professional in life in general. I also want to include Iiuri Karpenko for his help and advice in developing my programming infrastructure and Roland Katz for sharing part of his results from his previous contributions in the same field (quarkonia), which helped in the development of the model that I developed my thesis. I want to thank all the Subatech staff (my fellow Ph.D. students, secretaries, and other researchers) for making my stay in the laboratory a pleasant period of my life. In particular, I want to thank the Committee de Suivi de thèse members, Gines Martínez and Yacine Methar-Tani, for monitoring my work and offering beneficial advice

Introduction

In recent years, the research in the field of heavy-ion collisions (HIC) has been characterized by intense work in the experimental area, carried out by the collaborations associated with the LHC collider (ALICE, CMS, LHCb, and ATLAS), as well as the collaborations associated with the RHIC collider (PHENIX, STAR,. . . ). In the theoretical and phenomenological area, the discussions about the different models and approaches that attempt to explain the HIC process continue to be very debated, despite the extensive experimental data that the different collaborations have collected. In particular, the debates of the models referring to the production of the different particles is by far not over.

Despite the heated discussions, some consensus has been reached regarding which observables are more suitable for exploring the dynamics of HIC. Among the best observables, the production of particles stands out as one of the most promising observables, especially those that originate from the initial phases of HIC. Therefore, creating theoretical and phenomenological models that are able to predict the production of particles in HIC is a research priority for the community of heavy ion collisions.

The main objective of this work is to present and discuss a phase space coalescence model that allows predicting and exploring heavy quarkonium (bound state of a heavy quark Q and its corresponding antiquark Q) production in HIC as well as in proton-proton collisions. Our particular interest in heavy quarkonia is based on the fact that an important fraction of quarkonia is produced in the early phase of a heavy-ion collision, then travels through the quark-gluon plasma (QGP), and finally leaves the plasma without reaching full thermalization (so they do not lose all their memory of the processes in which they participated). By studying quarkonium formation with our model, we can obtain information of how the different interactions of quarkonium systems (Q Q) with the particles in the medium modify the production probability. Including the effect of QGP variables such as temperature and macroscopic characteristics of the plasma such as collective effects (flow).

Contents

To accomplish that, the thesis will be organized in the following way : The first chapter presents the main theoretical concepts and general definitions related to heavy-ion and proton-proton collisions. It Includes a short characterization of the event generator EPOS (Energy conserving quantum mechanical approach based on Partons, partons ladders, string, Off-shell remnants, and Splitting of partons ladders), and the Monte Carlo EPOSHQ version in which our model was inserted, in particular, a brief description of the conception of its initial state and hydrodynamic phase.

In the second chapter, several models and experimental results obtained by the Quarkonia community will be reviewed. Many of the experimental results and the models presented serve as a background and references for our work, together with a short description of the most popular tools and methodologies employed in the QGP research area.

In Chapter 3, the Remler formalism will be introduced. This formalism is the basis of our entire model. Our modifications and implementation details of the Remler model to the HIC collision scenario will also be discussed. Starting with the applications of our model, in Chapter 4, the charmonium production in proton-proton collisions at LHC energies will be presented and discussed.

Finally, in chapter 5, the results obtained for the production of charmonium for AA (PbPb) collisions will be exposed. We will close this work with the conclusions and considerations derived from the main results and the implementation of the model in the LHC environment.

Chapter 1

Fundamentals of quarkonium production in proton-proton and heavy-ion collisions

One of the most exciting aspects of high-energy collisions (proton-proton and heavy-ion collision) is the production of quarks of the second and third families, especially the so-called heavy quarks (c, b and t). As the stable nuclear matter is made up of only the members of the first family (u, d), except for compact objects (like neutron stars), high-energy collision experiments are (for the moment) our only ticket to study the production of these exotic(heavy) quarks and of the hadrons that they can form.

Quarkonium Production in pp Collisions

Proton-proton (pp) collisions are a crucial prerequisite to understand and quantify the effects of the hot and dense medium created in nucleus-nucleus collisions, as well as the Cold Nuclear Matter (CNM) effects. In particular, the study of the production of particles in pp (see results in figure 1.1) and its comparison with the cases pA and AA can directly point out the main differences between the three reactions (pp, pA and AA). The hadronic production of a quarkonium state (the one produced through interactions of quarks and gluons) is already intrinsically complex [START_REF]Measurement of J/ψ production in pp collisions at LHC energies with ALICE[END_REF], due to the different energy scales of the processes that take place during the formation of the bound state.

Nevertheless, they can be organized in two stages (adopting a simplified picture):

Chapter 1. Fundamentals of quarkonium production in proton-proton and heavy-ion collisions First, the heavy quark pair (HQ-pair) is produced in the parton-parton elementary interactions (hard scattering process). This is followed by the bound state formation (quarkonium) from the HQ-pairs, a process that occurs at a lower energy scale (soft process). The hard processes, responsible for the HQ-pairs' primordial (initial) production, are the gluon fusion and the annihilation of light quarks-antiquarks at leading order (LO) (see Fig. 1.2). Among the most promising models that try to explain the production of quarkonium are the color singlet model (CSM), the color evaporation model (CEM), and the Non-Relativistic Quantum Chromodynamics (NRQCD). The (CSM) conceives the formation of the quarkonia states as the creation of on-shell HQ-pairs from hard scattering processes described by perturbative QCD 1.1. Quarkonium Production in pp Collisions 13 (pQCD) with the same quantum numbers as the quarkonium state towards which they "evolve", where only the singlet color states are taken into account.

As for CEM, it is a model based on the quark-hadron duality, postulating that the cross section of a given quarkonium state is proportional to the cross section of the constituents of the pair, regardless of energy, transverse momentum, and rapidity.

And finally, NRQCD or color octet model (COM); in this model contributions from octet states of color are also taken into account, unlike CSM. Transitions between octet states to singlets states are treated as non-perturbative process employing an expansion in relative velocity powers between the HQ-pairs.

Color Singlet Model

The CSM assumes that the quantum state of the pair does not evolve between its production and its hadronization, neither in spin nor in color [START_REF] Braaten | Production of Heavy Quarkonium in High Energy Colliders[END_REF][START_REF] Lansberg | New Observables in Inclusive Production of Quarkonia[END_REF]. Although a gluon needs to be emitted in order to form the bound state, gluon emission from the heavy quarks is suppressed by powers of the α s (m Q ) coupling constant. If one further assumes that the quarkonium is a non-relativistic bound states with a highly peaked function in the momentum space, it can be shown that the partonic cross-section for quarkonium production should be expressed as that for the production of HQ pair with zero relative velocity v, in a color singlet state [START_REF] Redlich | Evidence for charmonium generation at the phase boundary in ultra-relativistic nuclear collisions[END_REF] and in the same angular momentum and spin state as the quarkonium to be produced, multiplied by the square of the Schrödinger wave function describing pair production in position space, taken at origin.

In the case of hadron-production, the cross-section of (for example) J/ψ can be computed by accounting the parton i, j densities inside the colliding hadrons, f i,j (x), and then factorizing the densities with the partonic cross-section of ij interaction to form a c c pairs

dσ[Ψ + X] = Σ i,j dx i dx j f j (x i , µ F ) f j (x j , µ F )d σ(µ R , µ F ) i+j→(c c)+X |Φ(0)| 2 (1.1)
The terms |Φ(0)| 2 can be extracted from decay-width measurement [START_REF] Lansberg | New Observables in Inclusive Production of Quarkonia[END_REF][START_REF] Redlich | Evidence for charmonium generation at the phase boundary in ultra-relativistic nuclear collisions[END_REF]. The model then becomes able to make predictions if we can provide some usual unknown values of non-physical factorization (µ F ) and renormalization (µ R ) scales of the HQ mass entering the hard part.

In the CSM, the states of quarkonium are interpreted as bound states of Q and Q Chapter 1. Fundamentals of quarkonium production in proton-proton and heavy-ion collisions that interact [START_REF] Braaten | Color and Spin in Quarkonium Production[END_REF][START_REF] Fritzsch | Producing heavy quark flavors in hadronic collisions-A test of quantum chromodynamics[END_REF] through a confining potential in a similar ways as presented in the figure 1.3. Using the J/Ψ as an example, it is identified as a bound state of a c c pair in a singlet state of color S 3 1 , the quantum number of the angular momentum being super index, and the subscript reflects whether it is a singlet [START_REF]Measurement of J/ψ production in pp collisions at LHC energies with ALICE[END_REF] or octet [START_REF]Measurement of quarkonium production at forward rapidity in pp collisions at √ s = 7 TeV,ALICE Collaboration[END_REF]. The only role of the gluons in the CSM model is to generate the potential that unites the Q and Q. In addition, the CSM provides a method to calculate the inclusive quarkonium production rate; it also provides the decay rate towards hadrons, and/or leptons and photons. Taking the J/Ψ as an example again, the expression for the decay of said particle into light hadrons according to CSM is: Γ(J/ψ) = Γ(cc(S 3 1 ))

|R J/ψ (0)| 2 4π (1.2)
Where Γ is proportional to the annihilation rate threshold in the S 3 1 state. The factor |R J/ψ (0)| 2 is the square of the wave function evaluated at zero distance between the c and c , which gives the probability for a c c pair to be found in a J/ψ state before the dissociation process takes place. In the CSM, it is assumed that the HQ-pairs that were produced in a hard scattering process will form a quarkonium state only if the Q Q pair is created in a color singlet state and the rest of the quantum numbers match those associated with the bound state to be formed.

Continuing with the case of J/Ψ, the HQ-pairs produced must be in the S 3 1 state, the formula for the cross-section of inclusive production for J/Ψ in the CSM is given 1.1. Quarkonium Production in pp Collisions by: σ(J/ψ) = σ(c c(S 3 1 ))

|R J/Ψ (0)| 2 4π (1.3)
where σ is proportional to the production rate of cc at threshold in c c(S 3 1 ) state. Again the factor |R J/ψ (0)| 2 is the probability that a point-like c c pair forms a J/ψ bound state with S 3 1 quantum numbers, and is given by the square of the wave function at the origin. FIGURE 1.4 -ALICE collaboration inclusive Υ(1S) p T -differential cross section [START_REF]Measurement of quarkonium production at forward rapidity in pp collisions at √ s = 7 TeV,ALICE Collaboration[END_REF] data point (red dots) at forward rapidity compared to a CSM [START_REF]Quarkonium production in proton-proton collisions with ALICE at the LHC,ALICE Collaboration and Philippe Rosnet[END_REF] at different orders (LO, NLO and NNLO)

Both Γ and σ factors can be calculated using pQCD based on the mass m Q and the coupling constant α s (m Q ). The only phenomenological parameter is the singlet wave-function of each spin multiplet: R J/Ψ (0) for the case S-wave of J/ψ and R χ (0) for the case χ c which are characterized by P-waves. The CSM provides predictions for the cross-section dependence with the polarization of the states of Quarkonium. The model also predicts that the ratios σ(H)/σ(J/Ψ) of the cross-sections of the different hadrons (H) with respect to the J/ψ as a function of the different quantum numbers J PC should vary strongly from process to process due to the angular momentum selection rules.

However, these last predictions proved to be wrong since significant variations in the ratios ( σ(H) σ(J/ψ) ) values have not been observed experimentally, and there has been Chapter 1. Fundamentals of quarkonium production in proton-proton and heavy-ion collisions very little or no evidence of spin asymmetries in the effective cross-section. The color singlet model has a strong predictive power; it has only one non-perturbative parameter for each state to be fixed. Since in color singlet model a hard gluon in order to form a quarkonium, the predictions are not reliable for larger values of the quarkonia meson fraction energy 0.8 < z < 1.0 (z ≈ E ψ √ s ). However, it was observed in [START_REF]J/ψ production at √ s= 1.96 and 7 TeV: Color-Singlet Model[END_REF] that when NLO (α 4 s ) and NNLO (α 5 s ) corrections are made to CSM, then it describes the data rather well (see for instance Fig. 1.4).

Color Evaporation Model

The CEM is based on the assumption that the probability that the produced HQ-pair forms a bound state is entirely independent of the color and spin states of the pair. This approach stems from assuming that the emission or absorption of soft gluons destroys any correlation between the color and spin states of the HQ-pairs and the quantum numbers of the final bound state.

This model is in line with the principle of quark-hadron duality. As such, the production cross-section of quarkonium is expected to be directly connected to the cross-section of Q Q pair production in an invariant mass region where its hadronization into a quarkonium is possible, that is, between the kinematics threshold to produce a quark pair, 2m Q and that to create the lightest open-heavy-flavor hadron pair, 2m H [START_REF] Redlich | Evidence for charmonium generation at the phase boundary in ultra-relativistic nuclear collisions[END_REF][START_REF] Amundson | Colorless states in perturbative QCD: Charmonium and rapidity gaps[END_REF]:

σ (N)LO (Ψ) = f Ψ 2m H 2m Q dσ (N)LO Q Q dm Q Q dm QQ (1.4)
The cross-section to produce a given quarkonium state is then supposed to be obtained after multiplication by a phenomenological factor F ψ related to a process-independent probability that the pair eventually forms a bound state. One assumes that some non-perturbative gluon emissions occur once the Q Q pairs are produced, and the emission of the non-perturbative gluon and the quantum state after the hadronization are essentially uncorrelated.

From the reasonable assumption that a 1/9 fraction of the Q Q pairs are in singlet configuration, the pairs in the suitable kinematics region of hadronization towards a quarkonium state can be counted as f ψ = 1 9 × 2J ψ+1

∑ i (2J i +1) based on a statistical approximation. Where J ψ is the spin of the quarkonium (charmonium) ψ and the sum over i runs over all the quarkonium (charmonium) states below the open heavy-flavor threshold.

The CEM does not attempt to cover (predict) the cross-section of decay, so it is a less ambitious approach than the CSM. However, it turns out to be much cheaper at the moment of calculating quarkonium production. It is important to note that the selection of 4m 2 H as the maximum value of the integration is entirely arbitrary, so s max = 4m 2 H also turns out to be a phenomenological parameter. In addition to providing an expression for the inclusive cross-section of quarkonium production, this model also predicts the cross-section ratios of the different quarkonium states. The model shows that the ratios between the different cross-sections are independent of the different production processes:

σ(H) σ(J/Ψ) = f H f J/Ψ (1.5)
The main consequence of the model is that the cross-sections are independent of the produced quarkonium polarization. For example, the fraction of J/ψ that is transversely polarized is predicted to be 2/3 for any process. This model captures the general features of the data. However, a quantitative comparison requires a model which incorporates phase space analysis. Nevertheless, it is important to mention that the CEM model still helps to predict rates of inclusive quarkonium production at large p T , while some other approaches, like for example the NRQCD (see later on this chapter), do not agree well [START_REF] Bodwin | Comparison of the color-evaporation model and the NRQCD factorization approach in charmonium production[END_REF] with the experimental data.

Even so, it is important to note that the traditional approach described in [START_REF] Redlich | Evidence for charmonium generation at the phase boundary in ultra-relativistic nuclear collisions[END_REF][START_REF] Amundson | Colorless states in perturbative QCD: Charmonium and rapidity gaps[END_REF] is not the end of the road for this model. Since a recent revision of the CEM resulted in an extension of the original idea in [START_REF] Vogt | Quarkonium Production and Polarization in an Improved Color Evaporation Model[END_REF], this approach can still contribute to the community.

Non-Relativistic QCD

The NRQCD factorization formalism has for objective the description [START_REF] Bodwin | Comparison of the color-evaporation model and the NRQCD factorization approach in charmonium production[END_REF] of the production and decay of heavy quarkonia. In this formalism, both production and decay are factored into short-distance coefficients and matrix elements of long-distance (LDMEs) [START_REF] Xu | Relativistic correction to color Octet J/ψ production at hadron colliders[END_REF], as represented by Fig. 1.5. The short distance coefficients solve the creation and annihilation of the HQ-pairs by a perturbative expansion of the strong coupling constant α s . The LDMEs, which represent the evolution of a free HQ-pair Chapter 1. Fundamentals of quarkonium production in proton-proton and heavy-ion collisions towards a quarkonium state, can be scaled by the relative velocity v between the members of the pair. They can be calculated using Lattice QCD or extracted directly from the experimental data.

The most important energy scales for the structure and spectrum of the quarkonium states are m Q v and m Q v 2 , where v 1. The approximate relative velocities values for the charmonia are v 2 ∼ 0.2 -0.3 and v 2 ∼ 0.08 -0.1 for the bottomonia.

The color octet model (COM) postulates that the HQ-pairs should be in color singlet (CS) bound state at long distances but that there may be octet (CO) state color at short distances. The COM was applied to cancel infrared divergences in the decay widths of the P-wave, and D-wave [START_REF] Braaten | Color and Spin in Quarkonium Production[END_REF][START_REF] Nlo Production | [END_REF] of the quarkonium states. In the FIGURE 1.5 -COM representation of J/ψ production in a pp collision. non-relativistic limit, the intermediate states (in the sense of time ordering of the perturbation theory) that contain HQ-pairs are suppressed and integrated out. In this way, the result obtained is an effective non-local action (at a distance of 1/m Q ). The field associated with HQ can be divided into two fields of Ψ spinors for the quarks and χ for the antiquarks. The coefficients of the operators are then tuned to reproduce the QCD employing a comparison between the Green on-shell functions computed in the QCD theory and those calculated according to the NRQCD Lagrangian.

At the desired accuracy, the α s (µ) and v 2 values determine which operators are maintained in the NRQCD and up to which loop-order they should be compared. Following this construction, the Lagrangian of the NRQCD has the same infrared behavior as that of the QCD, which implies that the result of the matching is independent of how m Q v, m Q v 2 and Λ scales are related. However, the relationship between these scales does influence which NRQCD operators should be maintained to obtain the desired accuracy. The general form of the NRQCD Lagrangian is [START_REF] Braaten | Color and Spin in Quarkonium Production[END_REF][START_REF] Bodwin | Rigorous QCD analysis of inclusive annihilation and production of heavy quarkonium[END_REF]:

L NRQCD = L 2 + L 4 + L gluon + ... (1.6)
The contribution to effective Lagrangian related to HQ has the following [START_REF]Non-relativistic effective theory for quarkonium production in hadron collisions[END_REF] form:

L 2 = Ψ + [iD 0 + -→ D 2 2m Q ]Ψ + 1 8m 3 Q Ψ + -→ D 4 Ψ + c 1 2m)Q Ψ +-→ σ • g -→ B Ψ+ c 2 8m 2 Q Ψ + ( -→ D • g -→ E -g -→ D • -→ E )Ψ + c 3 8m 2 Q Ψ + ( -→ D × g -→ E -g -→ D × -→ E )Ψ +... + charge-conjugated-for-the-antiquarks, (1.7)
where the first term describes a non-relativistic field according to Schrödinger's theory represented by the kinetic energy (Ψ + -→ D 2 2m Q )Ψ and a fixed term expressed as Ψ + (iD 0 )Ψ. The rest of the four terms were included to reproduce the Green functions of on-shell particles in the QCD theory up to order v 2 , including the hard gluon radioactive corrections expressed in terms of the chromoelectric strength fields E and B. To reproduce Green's functions with 2n external HQ fields, NRQCD must contain local operators with 2n quarks fields. The generic form (in our case, it has four quark fields) will be:

L 4 = ∑ i d i m 2 Q (Ψ + k i χ)(χ + k i Ψ)
+ four-quarks-scattering-operators, (1.8) where k i (k i ) are matrices whose indices run over color and spin and may also contain factors related to the -D/m Q derivatives especially in the case of high-dimension operators. The annihilation coefficients d i are in general complex, as their imaginary parts describe the decay processes of the quarkonium states. The Lagrangian of the NRQCD L 2 + L gluon coincides with the Lagrangian of the effective theory of heavy quarks (HQET) [START_REF] Golden | On the Hilbert space of the heavy quark effective theory[END_REF][START_REF] Balk | Heavy Quark Effective Theory at Large Orders in 1/m[END_REF]. Nevertheless, the two theories are different because the power counting schemes are different. In the Lagrangian of the HQET, the ordering is strictly by dimensions, which leads to a suppression of the kinetic operator (Ψ + -D 2 /m Q Ψ) in Λ/m Q factors (stationary description of the HQ). While in the NRQCD it is true that D ∼ m Q v and D 0 ∼ m Q v 2 indicating that the kinetic term cannot be neglected. The NRQCD has proven to be one of the most successful models in providing predictions for quarkonium production in a collision (1.6). It counts on a deeper theoretical base, including a Lagrangian one that describes the bound state formation process.

Color Glass Condensate

Thanks to the results reported by the HERA collaboration, we know that the densities of partons are very high at high energies or what is equivalent to small-x [START_REF] Iancu | The Color Glass Condensate and High Energy Scattering in QCD[END_REF][START_REF] Weigert | Evolution at small x bj : The Color Glass Condensate[END_REF]. This allows a description of this regime through the use of weak coupling physics. Therefore, high-energy QCD can be described as a many-body theory of partons that are weakly coupled, although with significant non-perturbative characteristics due to a large number of partons (color glass condensate). In this regime, the momentum distribution of gluons is peaked at some characteristic momentum value (the saturation scale Q sat ). However, as the energy increases, gluons are forced to occupy higher momentum states (as depicted on Fig. 1.7, left panel), which causes that the coupling constant becomes even weaker. The density of gluons is saturated at low values of α s (1/α s 1), and even resembles a multiparticle Bose condensate state.

In this picture, gluons evolve in large time scale with an upper natural limit of 1/Q sat [START_REF] Hirano | Mass ordering of differential elliptic flow and its violation for phi mesons[END_REF]. The gluons with small-x are classical fields produced by the gluons with higher values of x, and on the other hand, the gluons with higher momentum have a FIGURE 1.7 -On the left side, the different kinetic regimes (represented as a function of Q 2 and x) of QCD evolution. On the right side, the gluon density function as a function of the transfer momentum k T ( [START_REF]Heavy Quark production in the Color Glass Condensate framework: proton-proton and proton-nucleus collision studies in the ALICE Dimuon Forward Spectrometer[END_REF]) .

larger lifetime than those gluons with natural momentum scale (1/Q sat ). This means that gluons with small-x can be approximated by a static field (glass) and also that different configurations of gluons can be treated as non-interacting ensembles of fields (property of spinless glasses). Under this condition, if one assumes that the rapidity distribution of the gluons in the hadron is approximately constant, then the distribution of gluons will have the form

d 2 N p T dp T dy ∼ 1 p 2 T (1.9)
For small p T ≤ Q sat values (as long as p T ≥ Q sat ), on the other hand, the produced partons see a coherent field produced by a distribution of sources which are neutral color on the scale 1/Q sat as indicated in figure 1.7 (where k T = p T of equation 1.9).

This induces a dependence on the distribution of gluons (assuming again that the distribution in rapidity is constant) of the form

d 2 N p T dp T dy ∼ 1 α s (1.10) conditioning that: dN dy ∼ π 2 Q 2 sat α s (1.11)
Chapter 1. Fundamentals of quarkonium production in proton-proton and heavy-ion collisions

Saturation and Color Glass Condensate

Gribov proposed a simple criterion for gluon saturation, Levin Ryskin in [START_REF]Heavy Quark production in the Color Glass Condensate framework: proton-proton and proton-nucleus collision studies in the ALICE Dimuon Forward Spectrometer[END_REF][START_REF] Albacete | Gluon saturation and initial conditions for relativistic heavy ion collisions[END_REF], which amounts to comparing the density of gluons in a nucleon or nucleus and the cross-section for the recombination of two gluons. The density of gluons per unit area of the nucleus [START_REF] Iancu | The Color Glass Condensate and High Energy Scattering in QCD[END_REF] (of atomic number A) is:

ρ A ∼ xG A (x, Q 2 ) A 2/3 (1.12)
where xG A (x, Q 2 ) is the gluon distribution of the nucleus. Gluon recombination, whose cross section can be estimated by

σ gg→g ∼ α s Q 2 (1.13)
becomes important when ρ A σ gg→g 1. Indeed, this condition means that the interaction probability of a given gluon with the other gluons of the nucleus becomes close to one. When this regime is reached, the recombination slows down the growth of the gluon density, leading to saturation. The previous condition can also be rewritten as Q 2 Q 2 s where Q s can be defined as

Q 2 s ∼ α s xG A (x, Q 2 ) A 2/3 . (1.14)
Partons with Q 2 > Q 2 s are not affected by saturation, and their distribution evolves in accord with linear differential equations such as Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) or Balitsky-Fadin-Kuraev-Lipatov (BKFL). Conversely, those with Q 2 < Q 2 s are in a saturated regime. Their distribution evolves non-linearly and must be handled by the evolution equation of Jalilian-Marian-Iancu-McLerran-Weigert-Leonidov-Kovner (JIMWLK). The schematic representation of the JIMWLK equations can be understood in a generic non-linear structure of the form

∂n(k T , Y) ∂Y ∼ c 1 n(k T , Y) -α s c 2 n 2 (k T , Y) (1.15)
where the first term on the right side of the equation accounts for the gluon splitting, and the second term stands for the recombination of two gluons into one. With Y = ln(1/x), and where the coefficients c 1,2 ∼ 1. One can see that the quadratic term (that comes from gluon recombination) is not important as long as the gluon density n is much smaller than 1/α s . When n reaches 1/α s , the quadratic term becomes important, and as it appears with a negative sign in the equation, it slows down the growth of n with Y.

The equation (1.15) must be supplemented by initial conditions given at some moderate value of x = x 0 usually chosen of the order of 10 -2 . For this kind of x 0 , a common choice of initial condition for a large nucleus is provided by the McLerran-Venugopalan (MV) model [START_REF]Heavy Quark production in the Color Glass Condensate framework: proton-proton and proton-nucleus collision studies in the ALICE Dimuon Forward Spectrometer[END_REF], which describes the nucleus as a large collection of independent randomly distributed partons. From the MV initial condition, the solution of the BK/JIMWLK evolution equation gives the nuclear parton distributions at any value of x smaller than x 0 .

The CGC framework proved useful in studying heavy ions collision at high energy because the gluon saturation physics controls the bulk of the particle production in these collisions. Moreover, it has been shown that inclusive observables can be factorized (at leading logarithmic accuracy so far) as a convolution of universal distributions of color sources representing the gluon content of two colliding protons/nuclei and an observable evaluated in the retarded classical field produced by these sources.

Heavy Ion Collisions

From lattice QCD calculations performed at vanishing baryonic density, we know that a transition from ordinary nuclear matter to a deconfined state should occur if temperature (and then energy density) is high enough. The so-called quark-gluon plasma (QGP) is created once the temperature rises above a critical temperature T c ≈ 150 MeV [START_REF]Heavy-Ion Collision-hot QCD in lab[END_REF][START_REF] Heinz | Evidence for a New State of Matter: An Assessment of the Results from the CERN Lead Beam Programme[END_REF] at zero baryo-chemical potential, or more generally when density of the system is larger than ≈ 0.5 GeV/fm 3 . There are many similarities between the phase diagram of a fluid (water) and the phase diagram of matter QCD (see Fig. 1.8). Similar to the water phase diagram, equivalent states can be obtained by compressing and heating the matter.

The different phases can be accessed experimentally in high-energy heavy-ion collisions (HIC), in which an strongly coupled and hot fluid of QCD matter with low baryonic density is created after the initial interactions. Experiments in the SPS and RHIC collaborations have shown that a QGP essentially behaves as a strongly coupled liquid. In [START_REF]Phenomenological Review on Quark-Gluon Plasma:Concepts vs. Observations,Roman Pasechnik and Michal Sumbera[END_REF], these observations were now reaffirmed at LHC energies. Although curiously, correlations were detected between the particles produced in smaller systems such as proton-lead collision in which the QGP is not produced. A collision between heavy ions at ultra-relativistic energies can be divided for study into different stages as a function of time (see Fig. 1.9; τ at the collision's center of mass ) frame. The first phase of the collision is called the initial stage and involves the nuclei's initial radiation, the initial elastic collisions between nucleons, including the hard and soft partons scattering interactions. This phase is defined by a high energy density in which the interacting regions of the projectile nuclei overlap. This stage would usually last until around τ ∼ 0.3 fm/c (at LHC energy scale).

After the initial phase of the collision, the second phase (quark-gluon plasma stage) starts due to the numerous collisions in the initial phase. If this happens, the nucleons can no longer sustain their bound states. The quarks and gluons, now free, interact strongly in a small volume (volume of the fireball). The intensity of the interactions leads to a rapid thermalization of the partons in the medium.

At around τ ≈ 10 fm/c (at LHC energies), the thermalized QGP freezes-out rapidly, the quarks and gluons entwine, forming objects of white color (hadrons), this is known as the chemical freeze-out phase. Once the relative abundance of the different particles created is established (when the inelastic collisions end), the elastic collisions still occur between formed hadrons. They continue to interact elastically until the end of the kinetic freeze-out phase happening at τ ≥ 10 fm/c.

If one wishes to study the QGP properties in detail, one should pay attention to phenomena that are not directly related to QGP (due to its short lifetime of the fireball) but come as a consequence of QGP phases. These experimental methods are known mainly in the community as probes.

Among the different probes, we have the light particles production and their p T spectrum. This observable contains information about the hadron gas stage at freeze-out. Furthermore, these relative abundances of the particles and their momentum distributions are affected by the thermal properties of the QGP through the collective phenomenon known as collective flow (more discussion on this regard will be provided later in this thesis).

Another observable are the correlations of the produced particles. Since during the QGP phase, the partons are deconfined and strongly interacting, correlations between the particles produced in the final state of HIC can be modified as compared to the pp case. These correlations are sensitive to many properties of the system, such as the specific viscosity (shear viscosity), the emitter source size (via quantum correlations of produced particles), and the entropy density ratio (via the azimuthal momentum anisotropy of produced particles).

From the electromagnetic sector, one can also extract important information from the QGP. In that regard, we can also count on the production of direct low-momentum photons and the production of di-leptons of low and intermediate masses that are sensitive to the internal dynamics of the QGP. In particular, di-leptons provide an experimental window to access thermal radiation information from QGP via the Q Q → γ * → e -e + process.

Chapter 1. Fundamentals of quarkonium production in proton-proton and heavy-ion collisions

Jet modification is another well-known probe. Jets are clusters of strongly correlated particles originating from very high momentum transfer interaction with high virtuality. Jets can be strongly modified due to their interaction with the medium (jet quenching) and are thus sensitive to the density and temperature of the medium.

Finally, we have the heavy quark sector, where the most significant promise for QGP exploration is located. First, we have the so-called open heavy flavor sector, which can offer clear insights into the medium's transport properties, such as the diffusion and longitudinal drag coefficients, while the formation of heavy quarkonium (bound states of heavy quarks with their corresponding anti-quarks) is also sensitive to the medium's temperature. It is also important to emphasize that through the study of the formation and comparison of the different quarkonium states' relative abundances, the inner physics (thermodynamics) of the QGP can be efficiently assessed.

Coalescence Models of Quarkonium Production

Coalescence-based models are a promising alternative to classical models 1 of quarkonium formation that will be discussed in more detail in chapter 3. The origins of the coalescence models go back to the studies carried out on the formation of clusters in nuclear reactions. Initially, the coalescence models were conceived only in the momentum space. It was based on the idea that the formation of a compound particle (e.g. (anti) deuterium) was only possible if the (anti) nucleons were contained in a sphere in momentum space of radius p 0 [START_REF]Review of the theoretical and experimental status of dark matter identification with cosmic-ray antideuterons[END_REF], that is, if the relative momentum between both nucleons is less than p 0 . So the coalescence condition could be summarized in the way that for each pair pn whose relative momentum ∆p = |p pr n | ≤ p 0 , will coalesce, forming the corresponding compound particle (deuterium). The coalescence momentum scale p 0 is a phenomenological parameter, determined through experimental data [START_REF] Gev | Production of Tritons, Deuterons, Nucleons, and Mesons by 30[END_REF]. Despite the simplistic vision of these initial models, if a certain analytical distribution in the momentum space was taken into consideration for the single production of nucleons, then a consequent analytical distribution could be obtained for the different composite particles that can be obtained through coalescence. For example, deuterium production can be are the mass, the kinetic energy and the differential yield per event of the composed particle d (in this case deuterium) .

From these primordial models, the need to work on phase space stands out since relatively small values of relative distance or momentum respectively do not constitute a sufficient condition to guarantee the formation of bound states. The quark coalescence model can successfully explain the abnormally large enhancement of baryons to mesons ratio at intermediate p T and the scaling of the elliptic flow of hadrons according to their valence quarks [START_REF]φ and Ω production in relativistic heavy-ion collisions in a dynamical quark coalescence model[END_REF]. The mathematical tool used to develop models based on phase space is the Wigner function defined as

W( x, p) = d 3 y -→ x + y 2 |Ψ Ψ| -→ x - y 2 e -i h -→ p • -→ y .
(1.17) Thus, the cluster probability formation is expressed as an overlap (projection) integral between the Wigner function corresponding to the cluster wave function and the Wigner function corresponding to the N-body distribution of those nucleons which coalesce into the cluster, often approximated as the product of N single-body distributions.

For the case of meson formation (that is, the formation of Q Q bound state), the probability distribution as a function of p T of the final particle (meson) can be written as follows [START_REF]φ and Ω production in relativistic heavy-ion collisions in a dynamical quark coalescence model[END_REF][START_REF] Rhic | Parton Coalescence and the Antiproton/Pion Anomaly at[END_REF][START_REF]J/ψ production from charm coalescence in relativistic heavy ion collisions,Bin Zhang[END_REF][START_REF] Greco | Quark coalescence for charmed mesons in ultrarelativistic heavy-ion collisions[END_REF]:

d 2 P Ψ dp 2 T = g Ψ 2 ∏ i=1 d 3 r i d 3 p i (2πh) 3 E i W q (r 1 , p 1 )W q (r 2 , p 2 ) × W Ψ (r 1 , p 1 ; r 2 , p 2 )δ(p T -p 1,T -p 2,T )
(1.18) where f Q, Q(x i , p i ) are invariant distribution function of the quark and anti-quarks, which depend on the space-time position and four-momentum, including spin and color degeneracy. The factor g Ψ needs to be included in order to take into account the probability of forming a colorless (white) meson. Also considering the spin combinations, we end up with g Ψ = 1 9 × 3 4 for the case of the J/Ψ. And finally, W Ψ stands for the Wigner function of the final compound object (meson).

Chapter 1. Fundamentals of quarkonium production in proton-proton and heavy-ion collisions

EPOSHQ Event Generator

The model presented and discussed in this thesis has been implemented and executed inside the calculation infrastructure of the EPOS(HQ) event generator. The EPOS part of the event generator aims to reproduce the production of light particles in proton-proton and heavy-ion collision. The production process is conceived in five stages or steps: The initial conditions are defined based on the parton-based Grivob-Regge theory (PBGRT) to describe the objects (Pomerons) created in this stage. So called pomerons are treated as parton ladders, with a certain saturation scale and described by the DGLAP equations. Following the initial conditions, the core-corona approach takes place. Under this approach, the fluid region (core) and the peripheral or hadron escape region (corona) are identified and separated based on the momentum and density of string segments.

The third stage is characterized by viscous hydrodynamics, while in the fourth stage, statistical hadronization-like dynamics describe the hadron formation at the end of the QGP (or mixed phase). Finally, the fifth and last stage consists of a hadronic cascade implemented from the URQMD scheme.

Given the interests of this thesis, I will only briefly discuss the theory behind the partonic model in EPOS, which is crucial (among other things) in the initial production of particles and the hydrodynamic phase, which significantly influences the evolution of the heavy quarks (HQ) once these are created. Further information and a more detailed description of these subjects, among others related to EPOS, can be found in [START_REF] Matiere | Formation of Quark-Gluon-Plasma: Understanding the energy and system size dependence, Gabriel Sophys[END_REF][START_REF]Event-by-Event Simulation of the Three-Dimensional Hydrodynamic Evolution from Flux Tube[END_REF][START_REF] Epos | Charm and prompt photon production with[END_REF]

Parton Model in EPOS

At a theoretical level, the partonic model and the hydrodynamic evolution of the QGP play key roles in the operation of EPOS. The standard parton model approach consist in formulating the cross-section of a hadron-hadron or nucleon-nucleon scattering as factorization of the partons distribution function f h1 and f h2 in the projectile and target respectively, and the parton-parton cross-section σ h1h2 . In this way, the inclusive cross-section for the production of partons or jets in the final state as a function of transverse momentum (p T ) above a certain Q 2 0 scale writes

σ h1,h2 incl = ∑ ij dp 2 T dx + dx -f i h1 (x +,p 2 T ) f j h2 (x -, p 2 T ) dσ ij dp 2 T (x + x -s) (1.19)

EPOSHQ Event Generator

29

Where x = Q 2 2Mν is the momentum fraction of the interacting parton, ν = E -E is the change in the energy of nucleon during the collision and M the mass of the nucleon. The factor dσ ij dp 2 T is the elementary cross-section of parton-parton interactions while i, j represents the flavors of the partons.

This simplified factorization is the result of the cancellation of the higher-order diagrams (multi-scattering diagrams). These cancellations are referred to as AGK cancellations ( Abramovsky-Gribov-Kancheli cancellations) [START_REF] Alexey | AGK in the parton model[END_REF], which inhibit the structures of multiple scatterings of a given parton during the reaction. The most noticeable manifestation of the absence of these multiple scatter structures is that at high energy, the inclusive cross-section in proton-(anti-)proton collisions, is smaller than the total cross-section.

The solution to this problem is usually the implementation of the eikonal method, capable of dealing with contributions of multi-scatterings. Using this method, the inelastic cross-section can be calculated [START_REF] Theory | [END_REF] from the inclusive cross-section of the form:

σ h1h2 inel (s) = d 2 b(1 -exp -A(b)σ h1h2 incl (s) m! exp -A(b)σ h1h2 incl (s) (1.20)
which represents the cross-section for m scatterings; A(b) is the proton-proton overlaps function. The disadvantage of this method for recovering the contributions of the multi-scattering structures is the absence of a clear path for implementing the process of a collision between nuclei (AB). Figure 1.10 shows a graphical representation of the generalization of the parton-based Gribov-Regge Theory (PBGRT) formalism. This formalism is based on Feynman diagrams divided into several constituents (partons). Each one carries a fraction of the momentum of the incident nucleon. The sum of all the momentum fractions is equal to 1; therefore, the conservation of the initial momentum is ensured.

It is also important to note that each nucleon separates into several interacting partons (participants) and a single remnant. Participants interact with the other nucleon's partons via elementary interactions (vertical lines in Figure 1.10), and the remnants are simply the part of the nucleon that does not interact. FIGURE 1.10 -Contribution to the elastic amplitude of a nucleus-nucleus collision, or more precisely two nucleons from projectile A integrating with two nucleons from target B, taking into account energy conservation [START_REF] Theory | [END_REF]. The energy of the incoming nucleons is shared between all the constituents

EPOS Hydrodynamics Basis

Another fundamental component of the EPOS code that is influential in our model is the collision's hydrodynamic stage. After a very fast thermalization (around 1 fm/c) and initial partonic interactions, quarks, and gluons, evolve, showing a high degree of collective and correlations with other medium particles. The EPOS event generator conceives this stage in a very realistic [START_REF]Event-by-Event Simulation of the Three-Dimensional Hydrodynamic Evolution from Flux Tube[END_REF] way by considering the following factors: → Initial condition obtained from flux tube approach compatible with the string model together with the consideration of the possibility of an initial non-zero transversal collective flow. → Event by event procedure, taking also into account the high irregularity of the spatial structure of the single events. → Core-Corona separation based on the consideration that not all the matter produced in the collision has already thermalized, supplemented by the use of an efficient code to solve the hydrodynamic equations in a 3 + 1 dimensional space. → Uses a realistic equation of state, capable of reproducing the lattice gauge results, including a cross-over stage between the QGP and the hadronic gas phases. → Complete consideration of the hadron resonance table guarantees the compatibility of EPOS results in this phase with the statistical models. Finally, after thermalization, the evolution of the system is described by the URQMD hadronic cascade.

Due to the flux tube initial conditions, the system expands very quickly. However, the realistic proposal of the equation of state used in the EPOS model guarantees that the flow (radial and elliptical) and the collective effects, in general, develop earlier compared to the case of the first-order approximation of the equation of state (EoS) [START_REF] Hirano | Dynamical modeling of relativistic heavy ion collisions[END_REF][START_REF] Aichelin | EPOSHQ a new approach to describe charmed mesons in pp,pA and AA collisions[END_REF]. Therefore, the rapid expansion motivated by the initial conditions does not affect the collectivity of the system.

The system later begins the hadronization process in the cross-over stage (around a temperature of T = 166 MeV); the beginning of this phase is also the signal of the end of the QGP or thermalization phase. The evolution of matter under hydrodynamic laws stops at this point. Given fixed initial conditions, the core evolves according to the equations of an ideal fluid. The conservation of energy and momentum can be expressed as:

∂ µ T µν = 0 (1.21)
where ν = 0, 1, 2, 3. While the conservation of net charges can be written in the form:

∂N ν k = 0 (1.22)
being k = B, S, Q for the baryonic, strangeness and electric charge respectively. For an ideal fluid, the energy and momentum tensor can be decomposed as follows:

T µν = ( + p)u µ u ν -pg µν . (1.23)
Meanwhile the net charges can be written as:

N µ k = n k u µ (1.24)
where u is the four-velocity of the local rest frame. The solution of the previous equations provides the profile of the space-time evolution of macroscopic quantities such as the energy density (x), collective flow velocity v(x), and the net densities of flavor n k (x).

Chapter 1. Fundamentals of quarkonium production in proton-proton and heavy-ion collisions Even so, an additional ingredient is needed to ensure the complete and exact solution of the system of hydrodynamic equations; in this case, it is the EoS. Employing the equation of state, we can have a direct relationship between the energy density and pressure p, and this equation must also meet the following requirements: → Flavor conservation through the conservation of the respective chemical potentials (µ B , µ S and µ Q ) → Compatible with lattice gauge results for the case of (µ

B = µ S = µ Q = 0)
The first elements to take into account for the construction of an equation of state with the previously exposed premises are the pressure of the gas resonance (p H ) and the pressure of an ideal plasma of quarks and gluons (p Q ). If we consider T c , the temperature at which (p Q ) and (p H ) intersect, then the resulting pressure can be written as:

p = p Q + λ(p H -p Q ) (1.25)
where the dependence of the λ parameter is given by the expression

λ = exp - T -T c δ Θ(T -T c ) + Θ(T c -T) (1.26) being δ 0 δ = δ 0 exp - µ B µ c (1 + T -T c 2T c ) (1.27)
From the pressure we can obtain the density of the entropy S as:

S = ∂p ∂T = S Q + λ(S H -S Q ) + ∂λ ∂T (p H -p Q ) (1.28)
as well as the flavor densities n i by means of the equation:

n i = ∂p ∂T = n i Q + λ(n i H -n i Q ) + ∂λ ∂µ i (p H -p Q ) (1.29)
Finally, the energy density is given for the expression

= TS + ∑ i µ i n i -p (1.30)
The equation of state selected for inclusion in the EPOS model is obtained for δ 0 = 0.15 (referenced in figure 1.11 as X3F), which is capable of reproducing the results of the lattice gauge for the case

µ B = µ S = µ Q = 0. FIGURE 1.
11 -Energy density versus temperature, for our equation-of-state X3F (full line), compared to lattice data [START_REF]Event-by-Event Simulation of the Three-Dimensional Hydrodynamic Evolution from Flux Tube[END_REF] (black points) and some first-order EoS approach ("Q1F").

The X3F tag refers to the cross-over stage, where 3 is the number of active flavors in the state equation. Figure 1.11 also shows the results corresponding to the first-order equation of state referenced as Q1F, which is not good for large temperatures.

EPOSHQ

The EPOS-HQ model is a natural step forward in the evolution of the EPOS event generator, which does not include the production and evolution of HQs. This framework was designed to explore the behavior of heavy flavors observables in pp, pA, and AA collision. This tool is the result of the integration of MC@HQ into the EPOS3 code ( [START_REF] Aichelin | EPOSHQ a new approach to describe charmed mesons in pp,pA and AA collisions[END_REF][START_REF] Steinheimer | First results of EPOS-HQ model for open heavy flavor production[END_REF]). In this iteration, the heavy quark production was incorporated into the initial state of EPOS, being created similarly to the light-quarks, from gluon radiation at the initial or final state, or coming from very energetic gluon splitting.

The production of the HQ also contains cold nuclear matter effects, modeled by a saturation picture. At the same time, the propagation of the heavy flavors through the medium is carried out employing an Ultra-Relativistic Quantum Mechanics Dynamics [START_REF]Azimuthal correlations of heavy quarks in Pb+Pb collisions at LHC[END_REF], allowing the HQ to suffer an elastic collision. EPOS-HQ is the final and more crucial piece of infrastructure in which our model's code is placed. It offers a very accurate and realistic HIC scenario in which the HQ will be produced and evolve, as we can see in the figure1.12 which testifies that HQ dynamics is bona fide under control in EPOS-HQ. Further details on this dynamics will be provided in section 2.2.3.

Chapter 2

Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions

A unique feature of QCD matter compared to other strongly coupled systems is that its fundamental theory specifies the interaction directly at the Lagrangian level. Thus we have a chance to understand how a strongly coupled fluid emerges from a microscopic theory that is precisely known. The high temperature achieved in nuclear collisions, combined with the precision achieved by the numerical solutions of lattice gauge theory, permits the calculation of equilibrium properties of the hot QCD matter without any model assumptions or approximations.

At the same time discovery of strongly coupled state of matter (QGP) [START_REF]Exploring the Phases of QCD at RHIC and the LHC,US Heavy-Ion Community[END_REF] poses many questions. How do its properties behave over a broad range of temperature and chemical potential? What is the smallest droplet of hot QCD matter whose behavior is liquid-like? What initial conditions lead to hydrodynamics behavior, and can they be extracted from the experimental data? Can the underlying degrees of freedom in the liquid be resolved with jets and HF probes?

Answering these and other questions will also depend on an intensive modeling and computational effort to simultaneously determine the set of key parameters needed for a multi-scale characterization of the QGP medium as well as the initial state from which it emerges. Chapter 2. Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions

How Can We Investigate the QGP ?

Generally, one can divide the properties of the QGP into bulk medium variables, bulk transport properties, electromagnetic response, and medium response to propagating partons. The study of each of those categories often involves theoretical modeling of the dynamic evolution of different probes and their interactions inside the dense matter, ending in comparison with the experimental data for extraction of the properties.

Jets as QGP Probe

.1 -Jet quenching in a schematic head-on nucleus-nucleus collision. Two quarks suffer a hard scattering: One goes out directly to the vacuum, radiates a few gluons, and hadronises; the other goes through the dense plasma formed [START_REF]High-p T Hadron Suppression and Jet Quenching[END_REF] in the collision, suffers energy loss due to medium-induced gluons radiation (characterised by transport coefficient q and gluon density dN g dy ), and finally fragments outside into a (quenched) jet When studying high energy collisions, one often must consider processes where quarks and gluons are produced in the final state. These final-state quarks and gluons can be produced in several processes at LHC energies, from QCD radiation of incoming partons to W, Z, and Higgs boson decays. However, these high energy quarks and gluons are not observed in the final state of the collision due to the hadronization process. Overall, the high energy partons produced by the collision appear in the final state as a collimated bunch of hadrons, and they can be seen as proxies to the high energy quarks and gluons produced in the collision.

One of the most interesting experimental signatures which are suggested for QGP studies was the suppression of high-transverse momentum particles, often referred to as "jet-quenching" (see Figure 2.1 for illustration), resulting from in-medium energy loss when hard scattering partons pass through the strongly interacting medium [START_REF]Jet Quenching at RHIC and the LHC,Yen-Jie Lee[END_REF]48,[START_REF]Review of jet reconstruction algorithms[END_REF]. The typical lifetime of these high-p T showers is of the order of t ∼ 1 Q ∼ 1 p T , where Q is the parton virtuality scale.

A jet is an object made of several partons, and for that reason, a reconstructive algorithm needs to be applied to analyze the jets. Furthermore, the reconstructive algorithm must be the same for a direct comparison between the theoretical and experimental results. Although, the reconstruction algorithms are very useful tools, it may also introduce some uncertainties [START_REF]Review of jet reconstruction algorithms[END_REF], from the possibility of the change in the leading particle (the parton with the highest p T ) to background fluctuations and detector efficiencies.

Open Heavy Flavor

Heavy quarks are considered as "heavy" for two reasons: First, in the particle physics context, their mass is larger than the typical non-perturbative scale of the QCD [START_REF]Open Heavy-Flavor Production in Heavy-Ion Collisions[END_REF], m Q Λ QCD , which enables the evolution of their production cross-section within pQCD. Second, in the context of QCD matter formed in ultra-relativistic heavy-ion collision (URHIC), their mass is larger than the typical temperature reached in the QGP environment, m Q T QGP . This implies that HQ productions are essentially restricted to the initial hard scatterings (with a relatively short formation time of t f ∼ 1 2m Q ), and their thermalization time (which is parametrically ∼ m Q T 2 ) becomes comparable to (or even larger) than the fireball lifetime.

A key feature of HF probes in URHIC is their wide range in terms of transverse momentum, enabling systemic investigations of several processes that occur in different p T regions with a single probe, as illustrated on Figure 2.2. At low p T m Q , the mass of the HQs' enables to treat their propagation in the QGP as "Brownian-like" motion, with relatively small momentum transfer of the order of the temperature q 2 ∼ T 2 [START_REF]Open Heavy-Flavor Production in Heavy-Ion Collisions[END_REF]. Since energy is suppressed by T/m Q , elastic interactions dominate, the HQ motion for small HQ energy because the ratio between radiation energy Chapter 2. Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions and elastic energy loss is of the order of ∆E rad ∆E elast ∼ γα S , being α S the strong coupling constant. This allow us to reliably describe the evolution of the HQ in the medium by stochastic Langevin process characterized by a long wave-length transport parameter.

On the other hand, at sufficiently high p T , the mass effect ceases, and HF observable degenerates with those for light flavors, dominated by radiated energy loss. However, at intermediate p T , the dead cone effect is expected to suppress the small-angle gluon radiation, suggesting an energy loss hierarchy of the type ∆E b ≤ ∆E c ≤ ∆E q ≤ E g , leading to less suppression of HF hadrons compared to light-flavors hadrons. Open HF probes offer a unique opportunity to systematically investigate the interplay of radiative and elastic energy loss mechanisms over a broad momentum region and identify the transition between the two. Even so, there are some challenges that must be addressed to optimize HF's capabilities as hard probes. Two of the main issues of this probe are that we are still missing the precise dependence of the HF diffusion coefficient in the QCD matter and how the non-perturbative interactions driving the collisional diffusion of HQ at low momentum give way to the perturbative regime of energy loss via gluon radiation.

Collective Variables

The large energy loss suffered by HQ in the QGP is an indication of their strong coupling to the medium which is dominated by light quarks and gluons. If heavy quarks interact strongly with the medium, heavy flavoured hadrons could inherit the medium flow. Flow signals the presence of multiple interactions between the constituents of the medium created in the collision. More interactions usually leads to larger magnitude of the flow and brings the system closer to thermalization [START_REF]Investigating formation and evolution of Z(3) walls and flow anisotropies in relativistic heavy-ion collision[END_REF]. The magnitude of the flow is therefore a detailed probe of the level of thermalization.

Flow measurement can also lead to an understanding of the initial geometry, the initial energy density profile and the dynamical properties, such as the viscosity to entropy density ratio (η/s) of the medium produced in HIC. However, sensitivity to both, initial geometry and (η/s), acts as double edge sword, and it is difficult to tightly constrain either of these using the traditional v n measurements. To address this issue, some new observables have been recently developed [START_REF]Experimental overview on flow observables in heavy ion collisions[END_REF], which are sensitive to one but not the other. In order to separate collective effects due to the η/s from the effects due to the ones of initial geometry, correlations between the different flows (v n ) prove to be good tools to disentangle both effects. Correlations between different orders of flows establish the presence of non-linear hydrodynamics response in HIC, where higher-order flow harmonics are derived from lower-order eccentricities tied to initial geometry fluctuations.

A few years ago, a significant improvement in the understanding of HIC physics took place when researches realized that the initial collision geometry is affected by large event-by-event fluctuations. Following that line of thinking, more sophisticated methods [START_REF]Elliptic Flow in Pb-Pb Collisions at √ s NN = 2.76 TeV, ALICE Collaboration[END_REF][START_REF]ψ elliptic flow measurement in Pb-Pb collisions at √ s NN = 2.76 TeV at forward, L. Massacrier[END_REF] to study and estimate the even-plane angles Ψ n , which were previously assumed to be constant for a given event, were developed, thus arriving to the conclusion that even-plane angles have an explicit dependence on p T and pseudo-rapidity η.

Electromagnetic Probes

Radiation of photons and dileptons has been proposed as one of the most efficient tools to characterize the initial state of HIC. Unlike hadrons, which are emitted from the freeze-out surface after undergoing intense re-scattering, photons come out from each, and every phase [START_REF] Feinberg | Quark-gluon plasma and hadronic production of leptons, photons and psions[END_REF] of the expanding fireball (see Figure 2.3 for typical mechanisms). Chapter 2. Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions

Being electromagnetic in nature, they do not interact strongly with the medium, and their mean free path is larger than the typical system size (∼ 10 fm). After their production, they do not suffer further interaction with the medium (α QED α s ) and carry undisturbed information about the circumstances of their production to the detector.

The photon produced in the QGP (real and virtuals) can be used to studies the temperature of the plasma. Among most important possibilities are the evolution of the system size by intensity interferometry [START_REF] Photons | Dynamics of Quark Gluon Plasma and Interference of Thermal[END_REF], momentum anisotropy of initial partons [START_REF] Chatterjee | Anisotropic flow of thermal photons at RHIC and LHC[END_REF], as well as formation time of the QGP using the elliptic flow of thermal photons [START_REF] Chatterjee | Elliptic flow of thermal photons and formation time of quark gluon plasma at RHIC[END_REF]. The sources of direct photons can be classified into different categories depending on their origin from the different stages of the expanding fireball formed after the collision; there are: (1) prompt photons, which originate from initial hard scattering, (2) pre-equilibrium photons, produced before the medium gets thermalized, (3) thermal photons from QGP, as well as by hadrons reactions in the hadronic phase. Although it is not possible experimentally to distinguish between the different sources, a theoretical estimation can be used to identify these sources of direct photons and their relative importance in the spectrum [59].

Quarkonia

The quarkonia production (in particular quarkonia suppression) have been pointed out historically as a signature of the existence of QGP. This probe offers some advantages, first, the HQ which creates the quarkonia state are produced at an early stage of the QGP evolution and then propagate through the medium, participating in numerous interactions and escape the fireball before reaching a full thermalization (some memory will be preserved, especially those coming from primordial production and manage to escape the medium before getting melted). Secondly, due to their large masses, one can employ a non-relativistic treatment at the centre of mass frame of reference of the quarkonia state. The later is a consequence of the small relative velocity of the quarks in the quarkonium states, providing a hierarchy of the form [START_REF] Rangel | Quarkonium production in 2.76 TeV PbPb collisions in CMS[END_REF]:

m Q m Q v m Q v 2 .
The latter means that Effective Field Theories (EFT) can be applied. In particular, the J/ψ and especially the Υ (which posses a larger mass and whose production has an even smaller degree of thermalization) in HIC was suggested to be a promising probe of the QGP as a deconfined medium through screening processes of the two quarks of the quarkonium leading to their suppression.

Quarkonia Suppression

The particles selected as probes must not emerge fully thermalized from the medium; otherwise, they will lose the memory of the process they went through. Good candidates chosen as probes are heavy quarks (HQ) due to their large mass ( m c ≈ 1.3 GeV and m b ≈ 4.2) which "delay" their thermalization with respect to light quarks [START_REF] Isayev | Heavy Flavor Probes of Quark-Gluon Plasma[END_REF] and their thermal relaxation time is longer than the average lifetime of the QGP). In particular, the bound states formed by pairs of HQ and their corresponding antimatter quark (quarkonium) have a special interest because a comparison between their final distribution and corresponding initial distribution provides direct information of the properties of the QGP thanks to the screening effects on the Q Q interaction potential. The second consequence of its large mass scale is that its production is restricted to the earliest stages [START_REF]Heavy-quark production in heavy-ion collisions[END_REF] of the heavy-ion collision (HIC). Because only with hard scatterings (perturbative processes) can the heavy quark-antiquark pairs be created.

In this section, a compilation of the physical processes involved in the formation and evolution of quarkonium during the entire HIC process will be presented. Chapter 2. Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions

Cold Nuclear Matter Effects on Quarkonia Production

In general terms, cold nuclear matter (CNM) effects in p-A collisions can be quantified through the nuclear modification coefficient R pA

R pA = 1 A d 2 σ pA dydp T / d 2 σ pp dydp T (2.1)
Studies aimed to investigate these nuclear effects have been carried out since Matsui, and Satz [START_REF]ψ suppression by quark-gluon plasma formation[END_REF] first postulated the influence of these effects on the subsequent production of quarkonium, and more lately by the RHIC [START_REF]Cold Nuclear Matter Effects on J/ψ Yields as a Function of Rapidity and Nuclear Geometry in d+A Collisions at √ s NN =200 GeV[END_REF] and ALICE collaboration [START_REF]Rapidity and transverse-momentum dependence of the inclusive J/ψ nuclear modification factor in p-Pb collisions at √ s NN = 5.02 TeV,ALICE Collaboration[END_REF][START_REF] Arnarldi | Experimental overview on quarkonium production[END_REF].

From the interpretation of the experimental data obtained in the studies , it is concluded that there are other effects besides the color screening (Hot nuclear matter effect) that modify the formation and evolution of the quarkonium even before they are created. Those nuclear effects directly influencing quarkonium production are: nuclear absorption or normal suppression, the Cronin effect, and the so-called nuclear shadowing. Despite acting before the Q Q pairs formation, all the effects mentioned above, condition the scenario for the HQ pairs production. That conditioning leads to a discrepancy (R pA = 1) with respect to the prediction made in binary collisions between free nucleons, as illustrated in Figure 2.4. In the following sections, the main characteristics of the cited CNM and their influence on quarkonium production will be presented.

Nuclear Absorption or Normal Suppression

The idea of quarkonium absorption corresponds to the final-state interactions of the different states with the cold nuclear matter (nucleons). In order to calculate the effective cross-section of the charmonium produced in p-A collisions, Glauber formalism can be used [START_REF]Charmonium production and nuclear absorption in p -A at 450 GeV[END_REF][START_REF]Proton-Nucleus Cross Section at High Energies,Tadeusz Wibig, Dorota Sobczynska[END_REF], leading to

σ pA = σ 0 σ abs db[1 -(1 -T A (b)σ abs ) A ]. (2.2) Equation (2.
2) represents the break-up cross-section of any object that is traveling through nuclear matter, T A (b) is the nuclear thickness for an impact parameter b. In the formalism, σ 0 represents the elementary cross-section of nucleon-charmonium interactions. The simplest parametrization and one of the most used for calculations is:

σ pA = σ 0 Ae -σ abs ρL (2.3)
with ρL denoting the average area of nuclear matter traversed by the state (pre-formed) charmonium from its point of production. The value of ρL is obtained from the expansion of the integration of the expression (2.2) in powers of σ abs and remaining with the first term of the expansion.

ρL = A -1 2 d -→ b [T A ( -→ b )] 2 .
(2.4)

So the expression (2.3) is an approximation of the (2.2) valid for small values of σ abs .

There is still an additional way to account for the effects of nuclear absorption on quarkonium production, through the "α parametrization":

σ pA = σ 0 A α (2.5)
which will be equivalent to the previous approximation ( [START_REF]Charmonium production and nuclear absorption in p -A at 450 GeV[END_REF]) if α = 1σ abs ρL ln A , that also works for small σ abs values.

Cronin Effect

It was observed back in 1975 [START_REF]Production of hadrons at large transverse momentum at 200[END_REF] that the high-p T hadrons were not suppressed in proton-nucleus collisions, quite the opposite they have a reasonably large production rate. This effect was described by James Cronin, who demonstrates that bound Chapter 2. Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions nuclei cooperate to produce high-pt particles. A comprehensive description of this process was presented in [START_REF]Production of hadrons at large transverse momentum at 200[END_REF]. In this contribution, the influence of the center of mass energy on the Cronin effect is discussed. The mechanism of multiple interactions with nucleons, which is considered responsible for the broadening of the transverse momentum of the partons depends on the energy through the coherence length l c

l c = √ s m N k T (2.6)
At low energies, high-k T (transverse momentum) partons are produced off different nucleons incoherently, while at high energies, this process becomes coherent. For l c smaller than the typical separation distance between nucleons [START_REF] Cronin Effect | Hadron Production off Nuclei[END_REF][START_REF]Heavy-flavour and quarkonium production in the LHC era: from proton-proton to heavy-ion collisions[END_REF][START_REF] Johnson | Broadening of Transverse Momentum of Partons Propagating through a Medium[END_REF], the parton interacts incoherently (interacts individually with the nucleons), similar to the case of pp collisions. The projectile components interacts incoherently (nucleon by nucleon) with the nucleons of the target nucleus when the coherence lengths is shorter than the typical distance between nucleons. However, the QCD factorization is violated in these cases of multiple scattering [START_REF] Bodwin | Effects of initial-state QCD interactions in the Drell-Yan process[END_REF]. Therefore the broadening of the transverse momentum caused by the initial(final) interactions should not be moved towards a modification of the nPDF. In the opposite limit (for l c > R A ), the factorization is applicable, and the interaction with the nucleons is coherent in this case. From the system of reference of the incident parton(nucleons), the nucleus is a cloud of partons, with small x values, which means that it is no longer possible to discern among the different nucleons present in the target. The average value of the transverse momentum of the gluons increases since the gluon density is saturated at small values of k T . An approximate way to quantify the effect of multiple collisions is through the analysis of the modification of the kinematics of hard scattering a + b → c + d. For example, in p-A collisions, the distribution in the longitudinal fraction of the momentum x a is related to the equivalent in nucleon-nucleon collisions of the form:

[φ a (x a )] pA = [φ a ( x a 1 -a )] NN (2.7)
where a = ∆E a E a is the energy fraction lost of the incident parton (a) prior to hard scattering. In the previous relations, the modification suffered by the momentum fraction (x s ) due to energy losses [START_REF]High transverse momentum quarkonium production and dissociation in heavy ion collisions,Rishi Sharma and Ivan Vitev[END_REF] was taken into account. If the distribution function φ a (φ b ) has a Gaussian form then the variance of the transverse momentum of the partons acquires an additive property

k 2 Ta pA = k 2 T NN + k 2 Ta IS k 2 Ta IS = 2µ 2 L λ a (2.8)
in which we can appreciate the contribution to the transverse momentum of the collisions with the nucleons. In the expression (2.8), µ is the typical transverse momentum transferred in parton-nucleon collisions and λ a is the mean free path of the incident parton a.

Nuclear Shadowing

The measurements of the structure function F A 2 (x, Q 2 ) in the deep inelastic scattering (DIS) between leptons-nuclei [START_REF] Liuti | Consequences of nuclear shadowing for heavy quarkonium production in hadron-nucleus interactions[END_REF][START_REF] Eskola | Scale evolution of nuclear parton distributions[END_REF][START_REF]Nuclear shadowing (Topical Review),Nestor Armesto[END_REF][START_REF]Heavy Flavor Dynamics in Relativistic Heavy-ion Collisions[END_REF][START_REF]Cold nuclear matter effects on J/ψ and υ production at energies available at the CERN Large Hadron Collider (LHC)[END_REF][START_REF] Piller | Phenomenology of Nuclear Shadowing in Deep-Inelastic Scattering[END_REF] indicate that the parton density functions of the protons inside the nucleus are different from the corresponding functions of in free protons

f i/A (x, Q 2 ) = A f i/p (x, Q 2 )
. Nuclear shadowing effects are often categorized according to the value of the ratio

R A F 2 = F A 2 AF p 2
(2.9)

where F A 2 and F p 2 are the nuclear and nucleon structure functions, respectively. According to the range of the longitudinal momentum-fraction x and R A F 2 , different process [START_REF]Charmonium production and nuclear absorption in p -A at 450 GeV[END_REF][START_REF]Heavy Flavor Dynamics in Relativistic Heavy-ion Collisions[END_REF] dominate the interactions between the incoming partons and the Chapter 2. Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions nucleus. That is, if (R A F 2 < 1 and x ≤ 0.1) we will be in the presence of nuclear shadowing, if (R A F 2 ≥ 1 and 0.1 ≤ x ≤ 0.3) the anti-shadowing takes place. On the other hand, for the large values of 0.3 ≤ x ≤ 0.7 and R A F 2 < 1 the EMC effect is dominant and finally for x → 1 and R A F 2 > 1 the effects are due to the Fermi motion of the nucleons inside the nucleus. The results in [START_REF]Heavy Flavor Dynamics in Relativistic Heavy-ion Collisions[END_REF][START_REF] Piller | Phenomenology of Nuclear Shadowing in Deep-Inelastic Scattering[END_REF] show that the nuclear shadowing i) increases with the decrease of x, a behavior compatible with the CGC picture of the gluon saturation at small x discussed earlier, ii) also increases with the mass of the nucleus (A) and iii) decreases with increasing Q 2 . In the region of small-x, and relatively small values of Q 2 , the structure function F 2 (x, Q 2 ) can be estimated [START_REF] Piller | Phenomenology of Nuclear Shadowing in Deep-Inelastic Scattering[END_REF] approximately as:

F A 2 (x, Q 2 ) = Q 2 (1 -x) 4π 2 α CEM σ γ * -A (2.10)
In which γ * is a virtual photon that interacts with the nucleus. It is important to note that the cross-section γ * -A is mostly characterized by the interaction of the hadronic component of the virtual photon with the nucleus. The hadronic component of a photon is postulated in the vector dominance model (VDM), and this component increases with the photon energy. During its trajectory, the virtual photon can fluctuate in antiquark-quarks pairs (hadronic component) and can interact via strong interactions with the nucleus partons.

The importance of this phenomenon can be divided into two approaches: First, from a theoretical point of view, it allows us to test our understanding of the QCD in the high energy regime, while from the experimental point of view, one can also study the behavior of partonic densities within the nucleus for high energies.

Hot Nuclear Effects on Quarkonia Production

In the study of QGP produced in HIC, hard probes and electromagnetic processes can be used to explore many of the properties of the dense and hot medium. From properties of the space-time distribution of some observable to bulk properties and bulk coefficient, the signature of the strongly interacting medium can be extracted through the interpretation of the experimental data and the comparison with the theoretical prediction. Due to the presence of uncoupled quarks and gluons, the interaction potential of quarkonium states will be diminished (color screening effect). This effect is strongly dependent on the temperature of the surrounding medium of the quarkonium so that the temperature behavior can be studied through the bound states of the quarkonium. Other effects -to be discussed in detail in the next sections -also affect quarkonium production in HIC, as illustrated in Figure 2.7. Chapter 2. Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions

Debye Radius

The electrical potential generated by a test charge placed in a medium (electromagnetic plasma) can be shown [START_REF]Heavy quarkonia in quark gluon plasma as open quantum systems[END_REF][START_REF]Color Screening Melts Quarkonium,Agnes Mocsy and Peter Petreczky[END_REF] to be

V(r) = - Q 4π 0 r e -r r D , (2.11)
r D being the Debye radius. For distances r > r D , the potential drops exponentially, which indicates that an electric dipole whose distance is larger than r D will no longer be bound (dissociation); they will be seen as free charges. Something similar happens in the QGP with the color charges (although, in this case, there are three different "electric" charges). The QCD has a confining property that states that at large distances, the interaction potential grows. The potential in vacuum, has a linear term in addition to Coulomb one of the form

V(r) = - α r + σr (2.12)
where σ is the string tension, from which a more generalized parametrization can be postulated at finite temperature in a deconfined medium (QGP) by making both σ and α depend on the temperature of the medium:

V(r, T) = α r e -µ(T)r + σ µ(T)
(1e -µ(T)r ).

(2.13)

For temperatures above T c , the confining term of the potential becomes extremely small, so only the screened Coulomb potential term will remain. The magnitude µ(T) is the inverse of the radius of r D , also known as the Debye mass.

The screening mass has been calculated for a relativistic plasma in which was considered a number of colors N c and a number of flavors N f with a coupling constant g. It was shown in [START_REF]Heavy quarkonia in quark gluon plasma as open quantum systems[END_REF][START_REF] Burnier | A gauge invariant Debye mass and the complex heavy-quark potential[END_REF] that:

µ 2 = g 2 T 2 3 (N c + N f 2 ) (2.14)
and for the QED, the screening mass is:

µ 2 = 1 3 e 2 T 2 (2.15)
and finally for a two flavors QCD and N c =3

µ 2 = 4 3 g 2 T 2 (2.16)
Taking into account that the confining component is very small in a dense and hot environment, the Q Q interaction potential will take the form:

V(r, T) = - α r e -µ(T)r . (2.17)
This potential considers the medium's action on the pair, which can be assumed as the total effective action that the quarkonium perceives. More recent approaches based on lQCD calculations [START_REF] Mocsy | Quarkonia in the Quark Gluon Plasma[END_REF] nicely support this µ(T)/T scaling, as displayed in Figure 2.8. As illustrated in Figure 2.9, the effective potential in the medium starts to prevent the formation of bound states, as soon as the size of the mentioned states is larger than r D = 1 µ(T) . In comparison with this, in the vacuum (or cold nuclear matter), all are bound states have a size smaller than r D and thus exist, so that the missing of those states will lead to decreasing in the final production (suppression), being the later an indicator of the medium temperature effects .

Sequential Suppression of Charmonium

As discussed in the previous section, the effective potential of quarkonium in the medium depends on the temperature. Various quarkonium states are expected to be sequentially dissolved due to the color screening effect, depending on their binding energies and the temperature of the surrounding matter. As a result of that, quarkonium can play the role of a thermometer [START_REF] Mocsy | Quarkonia in the Quark Gluon Plasma[END_REF][START_REF]Quarkonium Binding and Dissociation:The Spectral Analysis of the QGP,Helmut Satz[END_REF] for the temperature reached in heavy-ion collisions and provides valuable information on the detailed properties of QCD in the high-density phase, as presented in Figure 2.10. The fundamental question is how to calculate the quarkonium dissociation temperatures. Once the

J/ψ(1S) χ c (1P) Ψ (2S) Υ(1S) χ b (1P) Υ(1S) T diss /T c 2.1 1.2 1.1 >3.0 1.2 1.2 ∆E[GeV]
0.64 0.20 0.05 1.10 0.67 0.54

Regeneration or Enhancement of Quarkonia

As presented in the previous sections (2.2.1 and 2.2.2), in the presence of a dense and hot medium (QGP), the quarkonium undergoes dissociation processes due to the action of the QGP. But in the same way, the dynamics of the QGP allows the existence of processes that contribute to the formation of quarkonium from the HQ present in the environment. The quarkonium regeneration has a quadratic dependence with respect to the initial density of c (b) produced through hard scatterings (gluon fusion and light quark-antiquark annihilation). The HQ in the medium can recombine due to thermal effects [START_REF] Zhou | Thermal charm and charmonium production in quark gluon plasma[END_REF] (thermal production). These regenerative processes are expected to play the leading role in quarkonium production at the FCC (Future Circular Collider) because they will be practically the only production mechanism as it is also expected that the temperature and density of the fireball generated at the FCC will annihilate the quarkonium initial production. Since the J/ψ meson was proposed [START_REF]ψ suppression by quark-gluon plasma formation[END_REF][START_REF]Heavy-flavour and quarkonium production in the LHC era: from proton-proton to heavy-ion collisions[END_REF] as a crucial observable to diagnose the production of QGP, it has been the inspiration of several theoretical models that seek to describe the evolution of J/ψ (υ) in the medium. One of the most notorious models designed for this purpose is the rate equation [START_REF]Heavy-quark production in heavy-ion collisions[END_REF][START_REF] Ganesh | pQCD approach to charmonium regeneration in QGP at the LHC[END_REF][START_REF] Rapp | Charmonium and bottomonium production in heavy-ion collisions[END_REF][START_REF] Rapp | Theoretical Overview on (Hidden) Charm in High-Energy Heavy-Ion Collisions[END_REF] model, which describes the temporal evolution of J/ψ abundances during the QGP phase. The corresponding equation for a system (fireball) of volume V FB where N J/ψ (τ) = n J/ψ (T)V FB (τ)n J/ψ (T) being the density of J/ψ at time τ -looks like:

dN J/ψ dτ = - 1 τ ψ [N ψ (τ) -N eq ψ (τ)] (2.18)
which is applicable as long as the J/ψ state exists and is well-defined, that condition will be fulfilled when τ -1 ψ = Γ ψ (τ) m ψ . For the right hand side of the equation, it is assumed that the surrounding partons (including c-quarks) are in thermal equilibrium. Chapter 2. Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions

The densities of the c-quarks and the J/ψ in equilibrium mesons are given by:

n c (T, γ c ) = n c (T, γ c ) = 6γ c d 3 q (2π) 3 f c (m c , T) n eq ψ (T, γ c ) = 3γ 2 c d 3 q (2π) 3 f ψ (m c , T) (2.19)
where f c/ψ (m c , T) are the thermal distributions of the partons/quarkonia in the fireball and γ c the fugacity. The other essential magnitude in the rate equation is the chemical relaxation time τ ψ . In the perturbative framework, the relaxation time can be determined by the cross-section of induced J/ψ breakup σ inel , which accounts for all the process that lead to the dissociation of the quarkonium bound state, via the convolution with the thermal distribution of the partons.

τ -1 ψ = d 3 k (2π) 3 Σ i=q,q,g f i (k, T)σ inel iψ .
(2.20)

Those rates [START_REF] Aronson | Collisional and thermal dissociation of J/ψ and υ states at the LHC[END_REF] are displayed as a function of time in figure 2.11, for a specific medium evolution. On the other hand, through a re-interpretation of the detailed balance equation for the J/ψ meson during its propagation in the medium, a variant of the rate equation [START_REF] Rapp | Charmonium and bottomonium production in heavy-ion collisions[END_REF][START_REF] Captain | Centrality and transverse momentum dependent suppression of υ(1S)[END_REF] can be extracted.

J/ψ + g ⇐⇒ c + c + X (2.21)
This variant illustrates more explicitly the contributions of the dissociative and regenerative processes separately.

dN ψ dτ = Γ rec N c N c [V FB (τ)] -1 -Γ diss N ψ (2.22)
Where Γ rec and Γ diss are the rates of recombination and dissociation, respectively. So if it is assumed that at the initial time of thermalization of the QGP (τ 0 ), the number of c (N c ) and of c (N c ) quarks are equal (N c = N c = N cc ), the analytical solution for the equation ( 22) in the case of N ψ < N cc will be:

N ψ (τ QGP , b, p T ) = (τ QGP , b, p T )[N ψ (τ 0 , b) + N 2 c c τ QGP τ 0 Γ rec (τ, b, p T ) V FB (τ, b) (τ QGP , b, p T ) dτ],
(2.23) Where N ψ (τ QGP , b, p T ) is the density of charmonium (bottomonium) states formed during the lifetime of the QGP, N ψ (τ 0 , b) is the number of J/ψ (Υ) produced initially at τ 0 . The magnitude (τ QGP , b, p T ) is the factor of suppression of the meson due to gluonic dissociation and collisional damping while the QGP still active. The suppression factor (τ) can be obtained as:

(τ QGP , b, p T ) = exp - τ QGP τ 0 Γ diss (τ , b, p T )dτ (2.24)
being τ 0 the time at which the J/ψ is created.

Heavy Quarks and Quarkonium Transport in QGP

During their propagation through the QGP, the HQ suffer frequent soft collisions with the medium partons. These soft interactions induce a Brownian-like motion for the HQ, which can be described by the Langevin transport equations [START_REF]Heavy Flavor Dynamics in Relativistic Heavy-ion Collisions[END_REF][START_REF] Li | Probing the transport properties of Quark-Gluon Plasma via heavy-flavor Boltzmann and Langevin dynamics[END_REF] (LTE).

dx i dt = p i E i dp i dt = F i drag + F i di f f (2.25) Chapter 2.

Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions

where

F i di f f = 1 √ dt C ij (t, p + ξdp, T)ρ j . The Gaussian noise ρ j follows a normal distribution P( -→ ρ ) = ( 1 2π ) 3/2 e -ρ 2 2 (2.26)
resulting in that: ρ i = 0 and ρ i ρ j = δ ij . Therefore, according to the previous properties, there is no correlation between the random forces between two different time scales

F i di f f (t)F j di f f (t ) = C ik C jk δ(t -t )
, which implies that the random kicks that HQ suffers are uncorrelated. For the other hand,

F i drag = -η D (p, T)p i , η D (p, T) being the drag coefficient. Parameter ξ ∈ [0, 1]
whose value depends on what momentum will be evaluated [START_REF] He | Relativistic Langevin dynamics[END_REF] in the argument of the covariant matrix C ij , e.g ξ = 0 for a pre-point (beginning of time interval), ξ = 1 2 mid-point or ξ = 1 post-point (end of time interval). The covariant matrix coefficients C ij can be represented as a function of the longitudinal (κ L ) and transverse (κ T ) diffusion coefficients, i.e

C ij (p, T) = √ κ L (p, T)p i p j + √ κ T (δ ij -p i p j ) (2.27)
Through equation (2.27), a relation between the drag coefficient η D and the coefficients κ L and κ T can be obtained [START_REF]Charmonium production and nuclear absorption in p -A at 450 GeV[END_REF][START_REF] Li | Probing the transport properties of Quark-Gluon Plasma via heavy-flavor Boltzmann and Langevin dynamics[END_REF], which is:

η D = κ L 2TE + (ξ -1) 1 2p ∂κ L ∂p + d -1 2p 2 [ξ( √ κ T + √ κ L ) 2 - (3ξ -1)κ T -(ξ + 1)κ T ], (2.28)
where d stands for the number of spatial dimensions that one may consider in the analysis, to obtain that expression (2.28), it was employed the rigorous definition of

η D . κ L and κ T η D = d p dt / p κ L = d (∆p z ) 2 dt κ L = d (∆p T ) 2 2dt
which describes the average momentum/energy loss and momentum fluctuations along the parallel (κ L ) and transverse (κ T ) directions to the direction of particle propagation.

It should be noted that if ξ = 1, a more compact expression can be obtained for In the results presented in the figure, the hard thermal loop (HTL) approximation was used to evaluated the soft interactions while the pQCD approach was employed to evaluates hard interactions with medium particles the relationship between the coefficients

η D = κ L 2TE - d -1 2p 2 ( √ κ T - √ κ L ) (2.29)
The sustainability of the transport model based on the Langevin equation is found in the high mass of HQ. This leads to assuming that interactions with the medium partons only produce small variations in the value of the momentum of the HQ (see results in figure 2.12). The evolution of the abundance and spectrum of the quarkonium states in the dense and hot medium can be tracked via the phase distribution f Q Q(r, q, t). The distribution function f Q Q satisfies the Boltzmann transport equation (BTE)

p µ ∂ µ f Ψ (p, t; p) = -ω Ψ Γ Ψ (p, t; -→ r ) f Ψ (p, t; p) + ω Ψ β Ψ (p, t; p) (2.30)
where Ψ = J/ψ, χ J/ψ , ..., Υ, ... stands for the possible quarkonium state presents in the thermal bath and ω Ψ = p 2 + m 2 Ψ is the total energy of the quarkonium. The coefficients Γ Ψ and β Ψ are the rates of dissociation and recombination respectively, for the other hand, the recombination rate accounting for the formation of the quarkonium. For a 2 → 2 process (via the inverse of the gluon-dissociation, Q + Q → Chapter 2. Characteristic and Main Features of Quarkonium Production in Heavy Ion Collisions g + Ψ), it takes the form:

β Ψ (p; r, τ) = 1 2ω Ψ d 3 k (2π) 3 2ω k d 3 p Q (2π) 3 2ω p Q d 3 p Q (2π) 3 2ω Q f Q (p Q ; r, τ) f Q( p Q, r, τ) × W gΨ Q Q(Θ[T diss -T(r, τ)](2π) 4 )δ 4 (p + k -p Q -p Q) (2.31)
The dissociation rate (Γ Ψ ) can be obtained by the cross section of the direct gluon-dissociation process

Γ Ψ = Σ i d 3 k (2π) 3 f i (ω k , T)v rel σ diss Ψ i (2.32)
Where v rel is the relative velocity with respect the incident particle (i.e gluon) and s = (p + k) 2 the center of mass energy of the collision. In order to secure the detailed balance, equation (2.32) must be equivalent to equation (2.20).

Low p T Quarkonia Production in Heavy Ion Collisions

Quarkonia states are produced before the beginning of the QGP's thermalization and are expected to partly dissociate into quark-antiquark pairs due to the screening of the strong interaction in the high color-charge density medium [START_REF] Michelangelo | Two Lectures on Heavy Quark Production in Hadronic Collisions[END_REF]. As mentioned before in this chapter, those dissociation processes can lead to a quarkonium suppression compared to the scaled (by the number of nucleon-nucleon collisions) production in pp. However, the suppression of the quarkonium states of charm and bottom quarks can also be due to CNM effects such as nuclear shadowing, energy loss, amount others effects. Even so, thanks to the relatively large yield of charm-anti charm pairs produced in HIC at high energies, charmonium states can also be created from the thermalized c and c via statistical production at the boundary phase [START_REF]Heavy-Ion Collision-hot QCD in lab[END_REF] of charms quarks in the plasma. The additional production mechanism is related to forming a hot medium and is particularly important for charmonium states given the higher number of the c c with respect to b b pairs. Recombination is expected to enhance the J/ψ yields, therefore mitigating or even counterbalancing the effect of the suppression mechanism.

Statistical Hadronization of Heavy Quarkonium

The statistical hadronization model (SHM) [START_REF] Andronic | Statistical hadronization of heavy quarks in ultra-relativistic nucleus-nucleus collisions[END_REF][START_REF] Braun-Munzinger | Evidence for charmonium generation at the phase boundary in ultra-relativistic nuclear collisions[END_REF] assumes that the charm quarks are produced in primary hard collisions and that their total number stays constant until hadronization. This approach assumes that the HQ are produce in primary hard collisions and that their total number stays constant until hadronization and that quarkonium mesons are then created during the hadronization process according to the available hadronic phase-space. Another important factor is the thermal equilibrium of the QGP at least near the critical temperature T c . This model was motivated by predicting the relative abundances of heavy hadrons produced in high energies HIC in terms of hadrons gas in chemical equilibrium. In the charm sector, the balance equation [START_REF] Braun-Munzinger | Evidence for charmonium generation at the phase boundary in ultra-relativistic nuclear collisions[END_REF] writes

N dir QQ = 1 2 g Q N th OQ I 1 (g Q N th OQ ) I 0 (g Q N th OQ ) + g 2 Q N th QQ (2.33)
where g Q is the fugacity factor, N Q Q is the number of initially produced Q Q pairs, and I n are the modified Bessel function. In a fireball of the volume V the total number of open (N th OQ = n th OQ V) and hidden (N Q Q = n th OQ V) charm hadrons is computed from their grand-canonical densities n th OQ and n th Q Q respectively. The balance equation (2.33) defines the fugacity parameter g c that accounts for deviations of HQ multiplicity from the value that is expected in complete chemical equilibrium. This can be understood in terms of the long time scales required to approach chemical equilibrium for heavy quarks, starting from a large number of charm (bottom) quarks produced via hard process during the initial stages of collision.

Considering the scenarios in which the formation of charmonium (bottomonium) is allowed to proceed through any combination consisting of one of the N c(b) quarks with one of the N c( b) antiquarks which results from the initial production of N Q Q pairs, one expects, for a given charm (bottom), that the probability P to form a charmonium (bottomonium) is just proportional to the number of available anti-charm (bottom) quarks relative to the number of light antiquarks

P c(b)→J/ψ(Υ) = N c(b) N u+d+s ≈ N cc(bb) N ch , (2.34)
where the number of produced charged hadrons is used to normalizes the number of light antiquarks. Since this probability is generally very small, one can simply multiply by the number of available charm quarks N c(b) to obtain the total number of J/ψ(Υ) expected in a given event.

In this framework, the absolute multiplicity of a quarkonium state (e.g. J/ψ see also results from the figure 2.13) can be estimated through the expression

J/ψ = g 2 c N J/ψ (2.35)
where N J/ψ is given by a thermal distribution of the form:

N j = d j V e P µ j /T 2π 2 Tm 2 j K 2 ( m j T ) ≈ d j Ve µ j /T ( m j T 2π ) 3/2 e - m j T (2.36)
where V and T are the volume and temperature of the hadron gas system. The variables m j ,d j denote particle masses, and degeneracy factors and K 2 is the modified Bessel function, and finally, µ j is the chemical potential of the particle j.

Prompt, Non-Prompt and Inclusive Quarkonium Production

As was discussed earlier in this work, the HQ production owns its origin to hard process in the early stages of the HIC. Even so, the formation of quarkonium states is controlled by non-perturbative processes (soft interaction) within a momentum transfer of the order of m Q v and their binding energy of the order of m Q v 2 , where v 2 ∼ 0.3 for charmonium and v 2 ∼ 0.1 for bottomonium. Those bounds states can be created with different quantum numbers (excited states) and each one of those states with their binding energy, radius, and lifetime.

Naturally, the ground states are the most stables among them, and for that reason, a significant rate of decays from the excitated states to the ground states should be expected. In that regard, a characterization of the different sources of quarkonium ground states (J/ψ and Υ) population need to be implemented to achieve a good understanding of the quarkonium production mechanism.

Charmonium Production

Inclusive J/ψ production is the sum of several contributions. In addition to the directly produced J/Ψ, the decays of the heavier charmonium states, such as the χ c and ψ 2S , also contribute to the observed yield, as illustrated in Figures 2.14 and 2.15. These two sources (direct and charmonium decays) [START_REF]Inclusive, prompt and non-prompt J/Ψ production at mid-rapidity in Pb-Pb collisions at √ s NN = 2.76 TeV,The ALICE collaboration[END_REF] are defined as prompt J/ψ. Since heavier charmonium states are less strongly bound than the J/ψ, they should be more easily dissolved in a deconfined medium [START_REF]Color Screening Melts Quarkonium,Agnes Mocsy and Peter Petreczky[END_REF].

On top of the prompt J/ψ production, there is an additional non-prompt contribution to the inclusive one at high center-of-mass energies collision, coming from the decay of the beauty hadrons. Only by subtracting the non-prompt contribution from the inclusive J/ψ yield, can one provide unbiased information on medium modification of the prompt charmonium.

Bottomonium Production

A significant difference in the bottomonium production with respect to the charmonium case, is that the former admits no non-prompt contribution from the feed-down of long-lived b-hadrons states, while the existence of the three sets of bottomonium states below the open-beauty threshold offers a wider variety of states that contribute to the Υ 1S state and can be studied. From the results of LHC, it can be seen that 70% of low-p T Υ are direct [START_REF]Heavy-flavour and quarkonium production in the LHC era: from proton-proton to heavy-ion collisions[END_REF], the second largest source comes from χ b . Meanwhile, at large p T (p T 20 GeV) less than half of the Υ are directly produced (see Figure 2.16 for more details). This can be explained by the fact that Υ(1S) state is the most tight bound state among all quarkonium states, implying that quite small suppression should be expected even at LHC temperatures. Meanwhile, bottomonium excited states are less tightly bound and, for this reason, more suppressed at low p T .

Low p T Inclusive Quarkonium Production

Thanks to Matteo Cacciari and collaborator's predictions (FONLL calculations), as well as the experimental results reported by the different collaborations ([104]), we know that quarkonium production is strongly favored at low p T in pp collisions as shown in the figure 2.17. This condition that the abundance in terms of p T is strongly concentrated in the region of low p T . The low p T quarkonium can be strongly affected by the medium because it has relatively less energy than those with high p T , which means that they will stay more time inside the fireball and the energy scale is closer to the one from the thermal gluons. The higher interaction with the medium can lead to an increase in the dissociation rate (2.20) so that almost all of the primordial quarkonium produced at low p T will be melted. According to the coalescence model's predictions in the low p T region, recombination is also more likely, especially in the final stage of the QGP when the temperature is not so high (thermal gluons do not have enough energy to destroy the bound states of quarkonium) and bound states are more easily formed. Therefore, despite losing most of the primordial quarkonium production, an enhancement of the production at low p T should be expected due to the recombination processes in the late QGP phase. The quarkonium produced through recombination processes, formed by the heavy quarks dissolved in the medium, will inherit parts of the medium's collective properties, such as correlations and flow. Although there are still discrepancies in how the quarkonia acquires these features from the medium, it is fair to say that recombination plays a crucial role in the generation of flow.

Collective Effects and Flow Coefficients

The flow at low p T signals the existence of multiple interactions between the constituents of the medium. The greater the number of interactions, the greater the flow measurement, and the system will move towards thermalization more quickly. Flow is usually described in the different theoretical and phenomenological models either as a consequence of the hydrodynamic behavior or as a microscopic transport of the QCD matter. In transport models, the flow depends on the medium opacity, either partonic or hadronic. In contrast, hydrodynamics is applicable when the mean free path of the particles is much smaller than the system dimensions. The most practical way to detect and study flow comes from observing the anisotropy in the momentum distribution of the particles produced with respect to the reaction plane. The reaction plane, being defined by the z-axis and the impact parameter vector, is illustrated by Figure 2.18.

The definition of the flow coefficients comes from the Fourier series expansion of the invariant distribution [START_REF] Snellings | Elliptic Flow: A Brief Review[END_REF]106] where E is the particle's energy, p is the three-momentum, p T is the transverse momentum, φ is the azimuth angle, y is the rapidity, and Ψ RP is the angle of the reaction plane. In the Fourier expansion of the distribution, the sinusoidal terms are all zero due to reflection symmetry with respect to the reaction plane.

E d 3 N d 3 p = 1 2π d 2 N p T dp T dy (1 + 2 ∞ ∑ n=1 v n cos(n(φ -Ψ RP ))), (2.37)
From equation (2.37), the flow coefficients have the form

v n (p T , y) = cos[n(φ -Ψ RP )] , (2.38)
where the dependence on p T and rapidity should be highlighted. The angular brackets indicate that the average is taken over all the particles, adding above all the events, in turn, must be taken in each bin (p T , y). The first flow coefficient (v 1 ) is called direct flow, which can be expressed as a function of the moments measured with respect to the reaction plane of the form v 1 = p x p T , while the second coefficient (v 2 ), called elliptic flow, can be computed as a function of the moments [107]

v 2 = p 2 x -p 2 y p 2 T .
(2.39)

The direct flow (v 1 ) is caused by the particles side ward motion within the reaction, while the elliptic flow (v 2 ) is due to the spatial geometry of the overlapping zone between the colliding nuclei. If the collision is not completely central,due to the differences in pressure gradients in the x and y directions, the spatial deformation will translate into a momentum distribution, so p x and p y are not generally equal, which implies that non-zero v 2 . The time scale for the conversion of the anisotropy in the pressure to that in momentum distribution measures the force of the medium interaction. A typical evolution is sketched in Figure 2.19. This anisotropy in the momentum distribution not only affects the light quarks (hadrons) but is also transmitted through interactions with bulk particles to the HQs, as shown in Figure 2.20. At the limit of multiple scattering where the momentum transfer in each interaction is small, the motion inside a thermalized medium can be treated as Brownian motion, which the Langevin equation can describe.

Low p T Quarkonia Production in Heavy Ion Collisions

65

This scenario can explain in a reasonable way how the HQs (and finally their bound states) obtain some flow from the particles of the medium through elastic collisions. In particular, the HQs that form hadrons (D-mesons, B-mesons and quarkonium) through recombination processes, are those that exhibit large elliptic flow, as a consequence of a higher degree of thermalization.

On the theory side, it is crucial to realize that state-of-the-art transport approaches like the one of TAMU, which are in excellent agreement with the nuclear modification factor of J/ψ, fail to reproduce their v 2 for p T larger then 4 GeV/c. This observation constitutes one of the major motivation for our own theoretical developments presented in the following chapters.

Chapter 3 Remler Formalism: A Coalescence Model For Binding Particles

Since T. Matsui and H.Satz postulated the suppression of J/ψ in a confined and hot medium, many attempts have been made to explain successfully the creation and subsequent evolution of the quarkonium in the medium. For the present work, it is the main interest to present a new model that allows for studying the production of quarkonium in the medium through the temporal evolution of its probability of formation and decay during its journey inside the fireball. The model uses coalescence in phase space to describe the formation of quarkonium states from the combination of Q and Q whose positions in phase space are close to each other. The probability for the formation of a bound state between a Q and Q pair will be evaluated by using the Wigner density function associated with the density operator expressed in the relative coordinates ρ(,

) = | |Φ Φ| | of the Q Q bound state. W Φ (, ) = d 3 ye -ipy r + y 2 |Φ Φ|r - y 2 . (3.1)
The ideas of our model were postulated by E. A. Remler [109]. He solved the time evolution equations for the quantum density operator using Wigner densities. This allows for the calculation of the rate of composed particle production. This formalism has been later successfully applied to the problem of deuterium formation in nuclear reactions [110]. More recently this method (in particular the phase space coalescence) has been applied with success to the production of charmonium in proton-proton collisions [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF]. Our main objective is to construct a model, based on the Remler formalism, that is applicable to the conditions created in HIC.

Fundamentals of Remler's Formalism

All information about a system can be extracted from the N-body density operator of the system (ρ N (t)). The density operator of the complete system obeys the von Neumann equation [110]:

∂ρ N /∂t = -i[H, ρ N ] (3.2)
where H is the Hamiltonian of the full system and can be decomposed as follows:

H = Σ i K i + Σ i>j V ij (3.3)
where K i is the kinetic energy operator of the particle i and V ij is the interaction between the particles i and j. Mesons, like a J/ψ, are two body objects described by the two-body density operator ρ J/ψ = |Φ J/ψ >< Φ J/ψ |. Φ J/ψ is the momentum eigenstate of a J/ψ formed from a Q and a Q. Thus

P J/ψ (t) = Tr[ρ J/ψ ρ N (t)], (3.4)
where the trace is taken over all N-body coordinates (which includes the Q and Q degrees of freedom), measures the probability of finding the Q and the Q in a subspace defined by ρ J/ψ at time t. It is difficult to use eq. (A.4) directly in some numerical implementations of the "cascade" type. Indeed, in this type of approach, the distance between the Q and Q quarks at the end of the QGP expansion is large with respect to the J/ψ radius and therefore P J/ψ (t) tends to zero. To overcome this problem we calculate the rate of the J/ψ formation and decay, Γ J/ψ (t), which is related to the probability to find a J/ψ at time T, P(T), by rewriting eq. (A.4) as

P J/ψ (T) = T 0 Γ J/ψ (t)dt (3.5)
with Γ J/ψ defined as

Γ J/ψ (t) = dP J/ψ dt = d dt Tr[ρ J/ψ ρ N (t)]
P J/ψ (T) = T 0 Tr[ρ J/ψ , ∂ρ N ∂t ]dt = -i T 0 Tr[ρ J/ψ [H, ρ N ]]dt,
where H is the Hamiltonian of the system. To obtain an useful expression for this rate we must first decompose the total Hamiltonian in the following way:

H = H 1,2 + H N-2 + U 1,2 (3.7)
where

H 1,2 = K 1 + K 2 + V 1,2 (3.8)
is the two particle Hamiltonian of Q and Q, H N-2 = Σ i K i + Σ j>i≥3 V ji is the Hamiltonian of the N-2 system and U 1,2 is the interaction of the particles 1 and 2 with the rest of the system.

U 1,2 = Σ i V 1,i + Σ i V 2,i . (3.9)
We perform this decomposition and replace the full Hamiltonian of the system in equation A. [START_REF] Lansberg | New Observables in Inclusive Production of Quarkonia[END_REF]. Profiting from the relations

[ρ J/ψ , H 1,2 ] = 0 (3.10) because ψ J/ψ is an eigenstate of H 1,2 and [ρ J/ψ , H N-2 ] = 0 (3.11)
because H N-2 does act only on the remaining N-2 particles due to the cyclic property of the trace we obtain

dP J/ψ (t) dt = Γ J/ψ (t) = -iTr[ρ J/ψ [U 1,2 , ρ N (t)]] (3.12)
Eq. (A.4) is the starting point of our approach. It allows to calculate the probability P(T) that a Q and a Q are in a bound J/ψ state at t = T by integrating the rate from t = 0 to t = T. P(T → ∞) is the probability that a J/ψ is observed in the final stage.

To make this possible we have to know ρ N (t). As a full quantum treatment of the evolution of ρ N is out of reach, one will resort to a semi-classical approximation The N-body Wigner density, W N , of the quarks in an expanding plasma is an object which cannot be calculated. In the past it turned out that many observable in heavy ion collisions can be well described if one replaces the N-body Wigner density by an average over classical phase space densities

W c W N ≈ W c N . (3.13)
Let us also here assume that the interaction between the N partons is of short range (as compared to the mean free path). This means that we consider that the QGP partons and the heavy quarks move on straight line trajectories between the collisions whose strength is given by cross sections. We can number the collisions by n and identify the scattering partners:i and j. We define as t ij (n) the time at which the n-th collision takes place which involves the partons i and j. This allows to calculate the momentum of particle i at time t as

p i (t) = p i (0) + Σ i =j Σ n Θ(t -t ij (n))∆p ij (n) (3.14)
where ∆p ij (n) is the momentum transfer in the n-th collision. ∆p ij (n) is equal -∆p ji (n). With this choice of W c N we can calculate the time evolution of the N-body density (eq.A.4)

∂W N (t)/∂t = Σ i v i • ∂ r W N (r, p, t) + Σ i≥j Σ n δ(t -t ij (n))× (W N (r, p, t + ) -W N (r, p, t -)) (3.15)
The first term arises from the straight line movement of the particles between the collisions while the second is due to the impulse received at the time t ij (n) when the n-th collision between particle i and j takes place. The δ(tt ij (n)) assures that a momentum transfer takes place exactly at the time of collisions. We can separate the change of W N due to collisions and straight line motion by writing

∂ρ N (t)/∂t = -iΣ j [K j , ρ N (t)] -iΣ j>k [V jk , ρ N (t)]. (3.16)
From the comparison between the equations (3.15) and (3.16) we find that

-iΣ j [K j , ρ N (t)] ≡ Σ i v i • ∂ r W N (r, r , t) (3.17)
and

-iΣ j>k [V jk , ρ N (t)] ≡ Σ i≥j Σ n δ(t -t ij (n))(W N (r, r , t + ) -W N (r, p, t -)) .
(3.18) Then substituting the equation (3.16 ) (the commutator between the potential and the full system density we end up

Γ Ψ (t) = ∑ i=1,2 ∑ j≥3 δ(t -t ij (n)) N ∏ k=1 d 3 r i d 3 p i W Ψ (r 1 , r 2 , p 1 , p 2 )× [W N ({r, p}; t + ) -W N ({r, p}; t -)] (3.19)
Hence, the elegant feature of Remler formalism is that only the scatterings of the Q and Q with the QGP medium need to be known in order to express the rate of quarkonium creation/destruction.

Non-Relativistic Wigner Functions

The Remler formalism presented in the previous section was the general frame to study cluster production in heavy ion collision around E kin = 1 AGeV as measured at the Bevalac at Berkeley. There it has been observed that in these reactions a large fraction of the nucleons are finally forming deuterons what has been not expected before.

This general approach has been transformed in numerical tool by Guylassy and Fraenkel who applied it so cascade calculations in which nucleons move on straight trajectories between collisions. For this they formulated the approach in terms of Wigner densities which are based on the same variables (r,p) as classical phase space densities. Later this approach has been been extended to Boltzmann-Uhling Uhlenbeck type of transport theories [122] and also for the formation of bound states of quarks and antiquarks (mesons) [START_REF] Showers | Jet Fragmentation via Recombination of Parton[END_REF].

The Wigner representation of any operator , e.g (O(t)), is related to its phase space Chapter 3. Remler Formalism: A Coalescence Model For Binding Particles coordinates representation by

W O ({r, p}; t) = dy 1 ...dy N (e -p 1 •r 1 ...e -p N •r N )× r 1 + y 1 2 , ..., r N + y N 2 |O(t)|r 1 - y 1 2 , ..., r N - y N 2 (3.20)
where r i and p i are the coordinates and momentum of the particles in the Wigner representation.

The N-body Wigner density, W N , of the quarks in an expanding plasma is an object which cannot be calculated. In the past it turned out that many observables in heavy ion collisions can be well described if one replaces the N-body Wigner density by an average over classical phase space densities

W c W N ≈ W c N . (3.21)
If the particles only scatter but do not have a potential interaction one talks about "cascade approach"; if the particles interact in addition via a mean field or by mutual two body interactions one talks about "Boltzmann Uhling Uhlenbeck approach" and "Quantum Molecular Dynamics", respectively. In the former two approaches the classical Wigner density is defined as

W c N = ∏ i h3 δ(r i -r i (t))δ(p i -p i (t)) (3.22)
as a product of point like particles with fixed momentum. In order to obtain physically meaningful results with classical Wigner densities we have to average over the time evolution of the classical Wigner densities. This is done by repeating the calculation many times. When it comes to the scattering of particles, the scattering angle is randomly chosen from dσ/dθ and then the collision is performed with that scattering angle. At the end of the calculations, when collisions have ceased, one averages over all finally observed particles obtained in these repetitions.

Two-Body Wigner Function

As the main interest of the model presented and discussed here is the study of heavy quarkonium production, a suitable Wigner function that can describe the two-body Q Q dynamics is needed. Due to their mass, most of the heavy quarkonium production is located in the low p T [START_REF] Rangel | Quarkonium production in 2.76 TeV PbPb collisions in CMS[END_REF] sector and thus a non-relativistic treatment is suitable for this region. Therefore under the considerations of long range forces are absent, non-relativistic dynamics (low p T sector) and taking into account that any interaction of either of the pair members with some external particle is treated as a punctual event in time (collisions), we may conclude that our heavy quarks behave like quasi-free particles that undergo a momentum change only when a collision occurs. This quasi-free evolution is only altered when a quark Q is close enough to one Q due to its potential Q -Q interaction. In this respect we have a situation which is close to that which Guylassy and Fraenkel assumed for the production of deuterons.

In this scenario one can assume that the wave function of heavy quarks in coordinate space wave function of the pair |Ψ i can be expressed as the product of a plane wave of the center of mass motion and the wave function of the relative motion |Φ i . The density operator of the HQ in the two-body state

ρ i = |Φ i Φ i | can be converted into Wigner density (3.20) W Φ i (r, p) = d 3 ye ip•y r - y 2 |Φ i Φ i |r + y 2 (3.23)
which allows to extend the wave function representation in coordinates space towards its dual space (the phase space), meaning that we can now represent each of the quantum state in a phase space representation through their the Wigner function associated with the density operator in coordinates space. The Wigner function has the normalization

∑ i d 3 rd 3 pW Φ i (r, p) = 1. (3.24)
From this we can derive the probability to find a

Q Q in a eigenstate i, if W 2 (r 1 , r 2 , p 1 , p 2 )
Chapter 3. Remler Formalism: A Coalescence Model For Binding Particles is the two-body Wigner density, by the expression

n i (R, P) = d 3 rd 3 pW Φ i (r, p)W 2 (r 1 , r 2 , p 1 , p 2 ) (3.25)
where the new coordinates R and P are related to the single particles coordinates by

R = r 1 + r 2 2 P = p 1 + p 2
representing the center of mass position and total momentum respectively whereas r and p are the relative coordinates associated with the pair measured in the center of mass.

r = r 1 -r 2 p = p 1 -p 2 2 .
(

If one assume that the probability density of each single pair is small, this approach can be extended to N-body system, by summing all of the possible Q Q pairs combination.

n i (R, P) = ∑ d 3 rd 3 p (2π) 3 W Φ i (r, p) ∏ j d 3 r j d 3 p j (2π) 3 W N (r 1 , p 1 , r 2 , p 2 , ...., r N , p N ) (3.27)
The product stands for all particles which are not part of the pair on which we are focus on at the moment. This means is j ≥ 3.

The initial relative wave function or Wigner density of the Q Q pair, can be calculated from first principles by solving the exact Schrödinger problem for a given potential. Thus taking inspiration from the harmonic oscillator approach, we can then approximate the Wigner function of two-body system to a double Gaussian-like function [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF] with respect to the relative distance r and the relative momentum p of the pair in the center of mass.

W Φ (r, p) = 8g (πh) 3 e -r 2 σ 2 -p 2 σ 2 (hc) 2 (3.28)
where g is a spin-color factors of a vector meson. The vector meson spin factor is 3.2. Non-Relativistic Wigner Functions 75 3 4 and the color factor for a color neutral Q Q pair is 1 9 , leading to a g = 3 4 × 1 9 . The factor 8 (πh) 3 is related to the normalization of the Wigner density (that was already established in equation (3.24)).

We have still to determine the Gaussian width σ of our two-body Wigner density W Φ . It is important to recall that the double Gaussian approximation (A.4) is not the exact solution of the two-body(meson) problem. But it can be calibrated from the exact Schrödinger equation solution a HQ pair in a Cornell potential. This exact solution [START_REF]RMS and charge radii in a potential model[END_REF] can then be used to compute the mean square radius

r 2 exact = Ψ(r)|r 2 |Ψ(r) d 3 r Ψ(r)|Ψ(r) d 3 r (3.29)
Then, deriving a relation between the mean square radius and the Gaussian width from the Gaussian approximation we find out that

r 2 Gaussian = r 2 e -r 2 σ 2 d 3 r e -r 2 σ 2 d 3 r = 3 2 σ 2 (3.30)
Combining now the equations (3.29) and (3.30) we will chose a value for σ in such a way that event with the double Gaussian approximation we will be able to reproduce the value of the mean square radius from the exact solution.

As our objective is to construct a coalescence model that operates in phase space, once the rate expression is obtained in any of the representations (positions or momentum), the next step is to move to the phase space using the associated Wigner function to the density matrix. The N-body Wigner density of the expanding medium will be replaced by the classical phase space density W c N , as discussed above. Now, expressing the probability rate in the phase space through the substitution of the probability densities by their corresponding Wigner functions we get

Γ 1,2 (t) = ∑ i=1,2 ∑ j≥3 ∑ n δ(t -t ij (n)) N ∏ k=1 d 3 p k d 3 r k h3 (W Φ (r 1 , p 1 ; r 2 , p 2 )× (W c N (r k , p k , ...; t + ) -W c N (r k , p k , ...; t -))) (3.31)
Paying more attention to phase integrals, one can realize that each factor W Φ (r 1 , p 1 ; r 2 , p 2 ) × W c N (r k , p k ; ...; t ± ) represents the projection of the two-body probability density over the full system probability density, which ultimately gives us the probability density Chapter 3. Remler Formalism: A Coalescence Model For Binding Particles of a bound state formation from particles 1 and 2 contained in the system of N particles.

Then defining the bound state formation probability at a given time t as

ω 1,2 (t) = N ∏ k=1 d 3 p k d 3 r k h3 W Φ (r 1 , p 1 ; r 2 , p 2 )W c N (r k , p k ; ...; t) (3.32)
which after substituting the semi-classical N-body Wigner density W c N (t) (which are delta function centered on the trajectories of the particles (r i (t), p i (t)) in the expression for the quarkonium probability ω(t) leads to

ω 1,2 (t) = W Φ (r 1 (t), p 1 (t); r 2 (t), p 2 (t)) (3.33)
meaning that the probability of formation of a bound state between particles 1 and 2 in a medium formed by N particles is equal to the two-body Wigner function where the coordinates and momentum are evaluated in the trajectory of particles 1 and 2.

Γ 1,2 (t) = ∑ i=1,2 ∑ j≥3 ∑ n δ(t -t ij (n))[ω 1,2 (t + ) -ω 1,2 (t -)] (3.34)
The ω 1,2 (t) which depends on the coordinates and momentum of particles 1 and 2, can be rewritten as a function of the position and relative momentum and the position and momentum of the center of mass.

ω 1,2 (t) = W Φ (r 1 (t), p 1 (t); r 2 (t), p 2 (t)) = δ(P -P(t))W Φ (r(t), p(t)) (3.35)
where r(t) = r 1 (t)r 2 (t) and p = 1 2 p 1 (t) -1 2 p 2 (t) are the relative position and the relative momentum, respectively, between the particles 1 and 2 evaluated in the classical trajectories and P(t) = p 1 (t) + p 2 (t) is the total momentum of the pair. It is important to recall that in equation (A.4) the delta function will be different from zero only in those cases in which one or both members of the pair (i = 1, 2) suffers a collision with any of the rest of the medium particles (j ≥ 3). Those interactions will of course modify the momentum of the particle and therefore the probability of forming a bound state. The latest equation summarizes the physical interpretation of the rate as the variation of the bound state formation probability when a collision happens (see figure A.1).

As was mentioned before, the rate Γ 1,2 gives us the probability rate by computing the difference in the probability of bound state (quarkonium) formation before and after any collision. But pushing forward our analysis one may argue that ω 1,2 (r(t -), p(t -)) represents a loss term, as it will always appears as with a negative sign in the rate equation, whereas ω 1,2 (r(t + ), p(t + )) can be identified as a gain term for the same reason (will appear always with positive sign). For each collision of a Q with a medium parton, both, gain and loss term, are calculated for all Q Q pairs in the system (and, vice versa, if a Q collides with a plasma parton we have to sum over all Q in the system). In the figure A.1, W -is the probability that a selected Q Q pair was a quarkonium state at tor in other words ((ω 1,2 (r(t -), p(t -))), shortly before the collision with a particle from the medium (light quark q or gluon g) . During the collision with a parton from the QGP Q gets a momentum transfer. W + is then the probability that the selected Q Q pair becomes a quarkonium state at t + that analogously we can identify as ((ω 1,2 (r(t + ), p(t + ))), shortly after the collision.

Going back to the relation between the rate and probability (A.4), we can evaluate the total probability that a HQ pair Q Q corresponding to particles {1, 2} forms a quarkonium state at any given time t by integrating the equation the rate (A.4) over Chapter 3. Remler Formalism: A Coalescence Model For Binding Particles

the time P 1,2 (t) = P prim (t 1,2 init ) + t t 1,2 init Γ 1,2 (t )dt (3. 36
)
P prim is the initial (primordial) probability of a quarkonium production of the Q Q, which also depends on the time (more discussion on this in chapter 5). The time integral of the rate accumulates all the variation in the probability of quarkonium formation of the pair Q Q during the time evolution due to any collision suffered either by Q or Q. The main characteristic of the rate defined in eq. (A.4) does not discriminate between the so called dissociation or recombination processes. This means that if during a collision some momentum ∆p Q, Q is exchanged between one of the members of the pair and a particle of the medium, this modification can lead to an increase or decrease in the probability of quarkonium formation. But in any case the real nature of the process will be reflected in our model by the sign of the rate.

Finally, the equation (3.36), which refers only to a single Q Q can be extended to a many Q Q pairs that are in the fireball at given time t. The later can be made by summing over all the possible pairs combination at a given time. The total probability of quarkonium formation at a given time can be calculated as

P tot (t) = N Q ∑ i=1 N Q ∑ j=1 P i,j (t) (3.37)
where N Q and N Q are the number of Q and Q that still inside the QGP as single particles (not being part of any bound state) respectively.

Relativistic Corrections to the Remler's Formalism

The relativistic formulation of quarkonium bound states is a difficult task, because we cannot anymore separate center of mass and relative motion. In addition, retardation effect can appear if the relative momentum of the pair is large in the center of mass and therefore and on top of that the times depend on the reference system.

Even so, in the appendices (A1, A2, and A3) we discuss a relativistic expression for the evaluation of the Wigner function in the laboratory reference frame (computational frame) as a function of the Wigner function evaluated in the center of mass. In these 3.3. Relativistic Corrections to the Remler's Formalism appendices we prove that the Wigner density in the laboratory frame can be evaluated as function of the Wigner function in the center of mass as show in the equation (A.39)

W i,u T,Φ ,y Φ (y, u T ; r, p) = 1 (2π) 3 δ(y -y Φ )δ 2 (u T,Φ -u T )W i,NR (r, p).
(3.38)

We so obtain the correct relativistic form of the Wigner function in any reference system. In eq. (3.38), the index i indicates a given quarkonium state Φ i , u T is the transverse component of the 4-velocity

u T,Φ = P T m Φ (3.39)
where P T is the transverse total momentum component of the Q Q center of mass, y and y Φ are the rapidity of the center of mass and the rapidity of the quarkonium Φ measured from the lab frame, and the index NR indicates that the Wigner function evaluated is non relativistic (A.4) and evaluated in the center of mass (taking advantage that due to the large mass of the HQ, generally the relative velocities in the center of mass are very small [START_REF] Rangel | Quarkonium production in 2.76 TeV PbPb collisions in CMS[END_REF]). Here is important to mention that only for those state in which we impose a well-defined center of mass 4-velocity and a well defined relative momentum with respect to the center of mass, we have been able to successfully derivative an prescription to evaluate our Wigner density at any system of reference as a function of the Wigner density in the center of mass frame. The later condition comes from the fact that even if we can always define a total momentum for the center of mass, the mass of the quarkonium state m Φ dependent on the relative momentum of the pair p as show in the equation (A.13) by the on-shell condition.

One can see from (3.38) that apart from the delta function which will tune the rapidity and the transverse component of the velocity with the corresponding values of the center of mass (rapidity and velocity), the real evaluation is carried out by the relative coordinates evaluated in the center of mass (cm) frame (r cm ,p cm ). Therefore, if we want to obtain the Wigner function in the computational frame (lab), we only need to boost the relative coordinates towards the center of mass (cm) and take advantage of the fact that in the center of mass we can use the non-relativistic expression of the function of Wigner W NR as long as the relative velocities of the HQ pairs are small.

Chapter 3. Remler Formalism: A Coalescence Model For Binding Particles

From that we can conclude that

W i,u lab T,Φ ,y Φ (y, u T ; r lab , p lab ) = 1 (2π) 3 δ(y -y Φ )δ 2 (u T,Φ -u T )× W i,NR (r cm (r lab ), p cm (p lab )) (3.40)

Appropriate Basis for the Quarkonium States in the QGP

As presented and discussed earlier in this chapter, Remler's original formalism was created to address the formation of nuclear light clusters in nuclear reactions from the resulting nucleons. As expected, the wave functions (probability density operators) that describe the particles involved in the process of formation of compound nuclei (particles) are the wave functions corresponding to the vacuum basis. That is, in equation (3.10) it is considered that the probability density operator ρ(r 1 , p 1 ; r 2 , p 2 ) corresponds to the two-body vacuum density operator. For that, the eigenstate does not change with time thus it will not contribute to the rate. This approximation is quite accurate for the nuclear reaction dynamics as well as pp collision in the QCD scenario. But is not good for the HIC situation, specially when HQ traverse the QGP.

Taking this into account, another important modification of our model with respect to the original formalism is thus the change of the states expressed on the vacuum basis to some local basis. The idea is to consider that the probability density functions that describe the formation of quarkonium states now also depend on the local temperature of the medium where the bound state is produced, implying that the local basis corresponds to the eigenstates of H 1,2 (T), where the temperature dependent potential is taken as the screened potential described in section 2.2.2 1 . This is the mandatory condition for (3.10) to hold in the QGP case. One can naively argue that for implementing such a change of vacuum to local basis, it is only enough to substitute in the equation (3.20) the two-body Wigner function, for another Wigner function which includes some dependency with the temperature. That is partially true because even if that rate depends on whether a collision between one of the members of the pair and a third particle happens or not, the probability density before and after the collision should consider also the medium effects. But that approximation will be also not good because as we are now considering the medium effects of the quarkonium bound state formation probability through the temperature, we also need to consider the change with time of the two-body bound state probability density ρ Φ . That is, go back to equation (A.4) and now consider that ρ Φ (r 1 , p 1 ; r 2 , p 2 , T(t)), which will lead us to an additional rate term (local rate term). Now, let us work out the expression for this local rate, as said before, we start from equation (A.4) and express the time derivative of the probability

Γ eff = Tr[ρ Φ (r, r , T(t)) ρN (t)] + Tr[ ρΦ (r, r , T(t))ρ N (t)] (3.41)
where the first terms is the "usual" Remler term, as we already know how that history ends, I will focus now only in the second term

Γ loc = Tr[ ρΦ (r, r , T(t))ρ N (r, r , ...; t)] (3.42)
Then after applying the trace integral, we obtain that the only integrals that are not equal to 1 are those related to r and r , so

Γ loc = d 3 rd 3 r ρΦ (r, r , T(t))ρ Q Q(r, r , t) (3.43)
where ρ QQ (r, r , t) = d 3 r d 3 r ...d 3 r n ρ N (r, r , r , ...; t) is the density operator of the N-body system after the trace operator was applied over the remaining N -2 particles which are not part of the pair, while ρ Φ is the density operator of the quarkonium bound state. Introducing the Wigner function through its relation with the density operator in the position space

ρ Φ (r, r ) = d 3 pe -i p•(r-r ) h W Φ (X, p) (3.44) being X = x+x 2 Γ loc = d 3 rd 3 r d 3 p e -i p •(r-r ) h W Q Q(X, p, t) d 3 pe -i p•(r r) h W Φ (X, r, T(t)) (3.45)
Then defining y = rr we can rewrite the Γ loc as

Γ loc = d 3 X 3 d 3 y d 3 p e -ip •y W Q Q(X, p , t) d 3 pe i p•y h ẆΦ (X, p, T(t)) (3.46)
Chapter 3. Remler Formalism: A Coalescence Model For Binding Particles which after work out the math (by integrating y ) leads to the final expression

Γ loc = (2πh) 3 d 3 Xd 3 pW Q Q(X, p, t) ẆΦ (X, p, T(t)) (3.47) Thus defining Γloc (X, p) = (2πh) 3 Ṫ(t)∂ T W Φ (X, p, T(t)) (3.48)
one get

Γ loc = d 3 Xd 3 pW Q Q(X, p, t) Γloc (X, p) (3.49)
Considering now the double Gaussian approach (A.4), which is the one that we use in our calculation

W Φ (X, -p) = W Φ (X, p) = 8g (πh) 3 e -( X 2 σ 2 + σ 2 p 2 h2 ) (3.50)
And remembering that the Wigner function W Q Q is nothing more that equation (A.4) applied only for two particles

W Q Q(X, p) = 2 ∏ i=1 h3 δ(X i -X i (t))δ(p i -p i (t))
is the semi-classical approach and the index i = 1, 2 will run just over the Q and the Q. Finally, including the temperature dependence of the Wigner function in the Gaussian width σ = σ(T), one gets

Γ loc = 8 Ṫ(t) dσ(T) dT ∂ σ e -( X 2 σ 2 + σ 2 p 2 h ) (3.51) = 16 (π) 3 σ(T(t))( X 2 σ 3 (T) - σ(T)p 2 h2 )e -( X 2 σ 2 + σ 2 p 2 h ) (3.52)
From the previous equation is now clear that in our double Gaussian approximation, the local rate term that emerges while considering a local basis (temperature dependent Wigner function), is due to the change in the Gaussian width with time. The physical interpretation of this is that the Gaussian width is directly related to the mean square radius of the bound state via

σ 2 (T) = 2 3 r 2 (T) (3.53)
3.4. Appropriate Basis for the Quarkonium States in the QGP

83

Thus one can argue that the broadening of the bound state dimension (mean square radius) with is responsible for that behavior. And that is the path through which temperature changes the bound state basis (local basis) in our approach. Furthermore, the time derivative of σ ( σ) is generally non-zero due to evolution of the medium temperature with time. During the QGP expansion the temperature changes dramatically, so for that some significant values of σ(T) are expected, especially during the early of the QGP evolution. The dependence between the Gaussian width σ and the local temperature T was obtained by fitting the results ([113, 114]) of the mean square radius of the ground state of charmonium (J/ψ) as a function of the temperature. The mean square radius ( r2 (T) ) was obtained from the solution of the Schrödinger equation with a Rothkopf potential. Figure A.2 illustrates the temperature dependence of σ with respect to temperature for the conditions of our simulation hereabove described. A corresponding value for the dissociation temperature T diss = 0.4 GeV is extracted accordingly, above which no J/ψ can be formed 2 .

Chapter 3. Remler Formalism: A Coalescence Model For Binding Particles

When working on the local basis the effective rate equation will be modified (A.4) (and therefore the probability equation as well (3.36) ) by the presence now of the local as indicated in the equation (A.62). The modify probability equation for a quarkonium state formation from a Q Q now has the form

P Q Q(t) = P prim (t Q, Q init ) + t t Q, Q init (Γ coll,Q Q(t ) + Γ loc,Q Q(t ))dt (3.54)
where Γ coll,Q Q is the collision rate contribution that can be calculated as was show in the equation (A.4), while Γ loc,Q Q is the rate contribution due to the change in time of the probability density (A.63) of the quarkonium state. Besides, When t → +∞, one has automatically T → 0 and the local basis thus correspond to the vacuum eigenstates which are measured experimentally.

Again, equation (A.4) can be extended to all pairs present in the medium at a given time by adding all possible combinations between Q quarks and Q antiquarks active (those inside the QGP, located in a region which the temperature is lower than the dissociation temperature and not being part of any bound state ) in the medium.

In-medium Q Q Interaction

Despite the Remler scheme was originally implemented in association with some cascade evolution of the N-body density taken at the semi-classical level, the genuine semi-classical evolution of W N turns out to have some role when particles are under the influence of the medium for times that are larger than the Heisenberg time (which is the typical time of Q and Q cycling each other in a given resonance). Indeed, as the rate equations ( (3.19) and (A.63) depend on the positions as well, it makes some difference whether Q Q stay close or not during the evolution time. Hence, Remler algorithm is just one half of the model, but it should be complemented by the other half which is a realistic semi-classical evolution of HQ in QGP medium including a state-of-the-art interaction between Q and Q.

The meta-potential U 1,2 in equation (3.7) represents the interaction of each particle that makes up the final composite state with the rest of the surrounding particles. This potential is what gives rise to the collision rate and is already well encoded inside MC@HQ generator. Even so, in our case, it is known that a quark-antiquark pair that is close to each other, interact via gluon interchange. This interaction -the potential V 1,2 in the internal Hamiltonian H 1,2 written at eq. (3.8) -is also present in conditions in the presence of an external environment (QGP), although with additional effects (screening effects). This interaction potential3 needs to be added to the semi-classical dynamics describing the HQ evolution. This was a considerable challenge in this thesis because the MC@HQ code had to be modified accordingly to run the HQ evolution under the designated interaction potential. Next subsections describe how this was achieved.

Relativistic Evolution of the Q Q Pairs

Each Q Q state will thus be propagated by solving the two-body problem in center of mass (cm) and then Lorentz boosting back and forth systematically from the cm to the computational frame (lab). This another difference between the current algorithm and the original Remler's formalism, which is evaluated in a unique frame (lab). The solution of the two-body problem was carried out in the framework of our relativistic correction to the original algorithm. This is because when we are dealing with a relativistic Q Q two major problems derived from the conservation of energy and momentum came to light. The non-conservation of energy during the evolution of a Q Q pair in a specific reference frame (eg the center of mass) can lead to tachyonic behaviors, which of course can destroy all computational effort. While on the other hand, for those pairs with very relativistic energies (momentum) can have also large relative momentum (even in the center of mass frame), and because of that some retarded interaction effects (Lienard-Wiechert like interactions) is expected.

We start out with the formulas of the relativistic kinematics which we will use to apply relativistic corrections to the non relativistic time evolution of a Q Q pair. The relativistic generalisation of the non relativistic equation of motion is given by

d dt (mv i ) = F i → d dτ (mu ν ) = dp ν dτ = K ν . (3.55)
where F i is the relativistic force and v i the non relativistic velocity, whereas u ν is the 4-velocity with

u ν = γ(v 1 , v 2 , v 3 , c) with γ = 1/ 1 -(v/c) 2 . K ν ,
u ν d dτ (mu ν ) = m 2 d dτ u ν u ν = d dτ mc 2 = 0 = u ν K ν (3.56) what yields K 0 = γ c F • v.
The Lorentz invariant proper time is given by dτ = dt/γ . Putting the components together we find

dp dt = v • F(r) F(r) (3.57)
The simplest relativistic modification to the movement of a Q Q pairs is to define a Lagrangian which contains the relativistic energy and the gluon fields in a Coulomb approximation

L = -γ -1 mc 2 + C r (3.58) with γ -1 = √ 1 -v 2 /c 2 and v 2 = ṙ2 + r 2 θ2
and m being the mass. With

∂L ∂v i = p i = γmv i (3.59)
and

H = p i v i -L = γmv i v i + mc 2 γ - C r = γmc 2 - C r = E (3.60)
we find H = E and hence a constant of motion. The equations of motion follow from the Euler -Lagrange equations

d dt ∂L ∂ qi - ∂L ∂q i = 0 (3.61)
for q i = {θ, r} and yield

d dt (mγr 2 θ) = 0mγr + m γ ṙ + C r 2 = mγr θ2 (3.62)
The first equation tells us that mγr 2 θ is a conserved quantity which can be considered as the relativistic analog of the non relativistic angular momentum and therefore

L = mγr 2 θ. ( 3
γr = -γ ṙ - L 2 mγr 2 d 2 dθ 2 1 r (3.66)
and the equation of motion is reduced to a trajectory equation

L 2 d 2 dθ 2 1 r -γ C m + L 2 r = 0 (3.67)
The last equation is a simple differential equation which can be solved for a given initial condition. If one wants to have ṙ one calculates r(θ) and uses later equation (3.65). After these relativistic considerations, anomalies related to relativistic behavior (tachyonic effects) have been cured. Throughout the previous analysis, an interaction potential V r Q has been taken into account because for the production of the bound state in hot and dense medium, some binding potential is need it. Although during the analysis for simplicity a Coulomb-like potential was assuming, the consideration of a different radial potential (Rothkopf-like potential) will only modify the form of the equation of motion (A.4).

Once we define the dynamics in the center of mass frame, we can use the strategy of running the HQ pairs evolution by systematically going to the cm, solving the dynamics as previously stated, updating the coordinates in the cm frame, then obtaining the updated coordinates in the lab frame to finally go back to the cm frame to repeat the process in the next time step. The drawback of such "back and forth" method is that small shifts in space-time coordinates are introduced at each step due to Lorentz transformations 4 . As those shifts are proportional to the Q Q 4-velocity, this approach is not suited for high p T pairs treatment. An alternative method would consist in implementing retarded Lienard-Wiechert like interactions. It was extensively investigated as well during this stage of the thesis. Although a Chapter 3. Remler Formalism: A Coalescence Model For Binding Particles very promising approach to this scenario, it comes with non conservation of the Q Q invariant mass due to associated classical radiation which is known to be inadequate at the quantum level. For that reason, we have decided not to follow this track.

Q Q Interaction Potential Effects

Having in our hands the tools to carry out the in-medium evolution of the HQs, the next step to do is to explore the effects of the Q Q interaction potential in the evolution of the pairs. The most direct way to study the effects of the interaction potential is to pay attention to the behavior of the relative coordinates (in particular to the average relative distance of the Q Q pairs).

The following figure (3.3) shows a comparison between three different scenarios: no Q Q interaction potential (blue dots), screened Q Q potential (green dots), and vacuum Q Q potential (red dots). In this study, we are interested in knowing how the instantaneous number of Q Q close pairs (number of Q Q pairs in a given time the (invariant) relative distance is lower than 1 fm) In Figure 3.3, we see that even with screening effects the Q Q interaction potential plays a significant role in the probability for quarkonium production, especially 3.5. In-medium Q Q Interaction 89 at the final stage of the evolution, by substantially increasing the number of Q Q pairs that stay close (with small relative distance). If there is no interaction potential (blue dots), we see that the number of Q Q that stay close diminishes quickly with time. This is quite important because in a fast expanding medium, the HQ tends to move away from each other, but adding an interaction potential help (especially the latest stage of QGP evolution) to the recombination (regeneration) process to take place and recover some of the initial (primordial) production that was melted by the medium dissociative power which is significantly stronger at early stages of the fireball evolution. Figure 3.3 also inholds the case where the vacuum potential is considered in the semi-classical dynamics. Although such scenario where an unscreened potential would apply in the QGP is clearly unphysical, it allows to gauge the the importance of the systematic uncertainties that are still paving the estimation of the Q Q potential.

The influence of having this Q Q interaction potential included in the semi-classical evolution, both for the collision rate and for the local rate will be studied rigorously in the chapter 5 of this thesis.

Chapter 4

Charmonium Production in pp Collisions

In this chapter the application of our model to proton-proton collisions (pp) will be presented and discussed. The main objective of this chapter is to introduce an auxiliary magnitude (the conditional probability) that allows us to calculate the probability that a single uncorrelated HQ pair can form a quarkonium state as a function of its transversal velocity (momentum). From this auxiliary magnitude and also assuming a factorization between the transverse magnitudes (velocity or transverse momentum of the center of mass) and the rapidity, the total production of a quarkonium states in pp collisions will be obtained.

During this chapter we used the notation which is summarized in the table relative momentum q q lab , q cm relative position r r lab , r cm "Laboratory" refers to the position and momenta in the computational frame. This is in the simulation of heavy ion collisions the nucleus-nucleus center of mass frame. "cm" refers to position and momentum in the rest system of Q Q system, which my form a quarkonium.

We introduce the center of mass 4-vectors

Total momentum P = p Q + Q Q2 generalize velocity U = γγu
Center of mass quantities are characterized by capital letters, while for relative or single particle coordinates we use small letters. In case of need we work with two sub-indices, one indicating a component and the other indicating an state(particle), Chapter 4. Charmonium Production in pp Collisions separated by a coma, like for instance U Φ,T , which indicates the transverse component of the generalized velocity of the quarkonium state i, measured in the center of mass.

For the quarkonium production in pp collisions we start out be introducing the conditional probability for the formation a quarkonium state from a non-interacting uncorrelated (Q Q) pair. The total production will be obtained by multiplying the conditional probability with all possible Q Q pairs. By assigning to each pair the center of mass rapidity we can finally obtain the rapidity distribution of the Q Q pairs.

Conditional Probability of quarkonium production

Our approach is based on the Remler formalism. It has already been applied to the study of quarkonium (charmonium) production by Song et al. [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF]. We will revisit the pp collision scenario in order to have a baseline for the calculations for heavy ion collisions where several modifications will be applied: Relativistic corrections, temperature dependent Wigner functions and collisions with medium particles. The knowledge of the results for proton-proton collisions is also crucial for obtaining the nuclear modification factor (R AA ) which will be studied in the following chapters.

We first have to introduce the concept of conditional probability. The conditional probability is the probability that a quarkonium state Φ is produced in a pp collision, in which an uncorrelated single Q Q pair has been created, assuming that the interaction between Q and Q is negligible in comparison with the momentum of the quarkonium in the center of mass frame. Thus, Assuming that the distribution of Q Q produced in a pp event is normalized to unity

d 3 p lab Q d 3 p lab Q d 3 r lab d 9 N Q Q d 3 p lab Q d 3 p lab Q d 3 r lab = 1.
(4.1)

p lab Q,
Q are the 3D momentum vectors of Q and ( Q), measured from the lab frame and r lab is the relative distance between Q and Q measured from the lab. Let us furthermore consider that

d 9 N Q Q d 3 p lab Q d 3 p lab Q d 3 r lab = d 3 N Q Q d 3 r lab × d 6 N QQ d 3 p lab Q d 3 p lab Q (4.2)
4.1. Conditional Probability of quarkonium production

93 with d 6 N Q Q d 3 p lab Q d 3 p lab Q = d 3 N Q d 3 p lab Q × d 3 N Q d 3 p lab Q . (4.3)
Knowing that d 3 p lab Q d 3 p lab Q = d 3 P lab d 3 q lab , where P = p Q + p Q is the center of mass momentum and q = 1 2 (p Qp Q) is the relative momentum of the pair we can rewrite equation (4.1)

d 3 P lab d 3 q lab d 3 r lab d 9 N Q Q d 3 P lab d 3 q lab d 3 r lab = 1. (4.4)
Employing the relation between the longitudinal component of the total momentum and the center of mass rapidity y

d 3 P lab = d 2 P lab T dP lab L = d 2 P lab T E lab Q Qdy (4.5) with E lab Q Q = E lab Q + E lab
Q being the total energy of the pair and E lab Q and E lab Q the quark (Q) and anti-quark ( Q) energy, respectively. P lab L and P lab T present the longitudinal and transverse momentum, measured in the lab frame. Thus we rewrite for the normalization condition [START_REF]Relativistic Wigner distributions[END_REF]

d 2 P lab T dyd 3 q lab d 3 r lab E lab Q Q × d 9 N Q Q d 2 P lab T dP lab L d 3 q lab d 3 r lab = 1 (4.6)
On the other hand, the relationship between P T and U T can be established through the equation

P 2 T = U 2 T s
where s is the Mandelstam variable, defined in our case as s = (p lab Q + p lab Q ) 2 , being p Q and p Q the 4-momentum of the quark and anti-quark respectively. The Mandelstam variable (s) in the lab frame will be a function of U T for a fixed relative momentum (q lab ) and that relation can be derive as follows

s = (p lab Q + p lab Q ) 2 = 2m 2 Q + 2p lab Q • p lab Q = 2m 2 Q + 2(E lab Q E lab Q -p lab Q • p lab Q ) (4.7)
then, taking advance of the inner product invariant p lab Q • p lab Q for computing that product in the center of mass where

p cm Q = -p cm Q and E cm Q = E cm Q , we get p lab Q • p lab Q = (E lab Q E lab Q -p lab Q • p lab Q ) = m 2 Q + (p cm Q ) 2 + (p cm Q ) 2
s = 2(m 2 Q + (m 2 Q + 2(p cm Q ) 2 )) = 4(m 2 Q + (p cm Q ) 2 ) (4.8)
Finally, knowing that p Q = P 2 + q and that in turn P cm = 0, it will be obtained that

s = 2(m 2 Q + (m 2 Q + 2(p cm Q ) 2 )) = 4(m 2 Q + (q cm ) 2 ) (4.9)
Then, returning q cm → q lab to the lab reference frame, we will have to consider the components of the vector q lab parallel and perpendicular to the vector U T . Therefore, decomposing the relative momentum as ||q lab || 2 = (q lab L ) 2 + (q lab T ) 2 and taking into account that q lab L is perpendicular to the plane where U T is located, we can finally derive for the relation between s, U T and q lab , by projecting the component of q cm vector parallel to U T and then boosting that component by the corresponding Lorentz factor.

s(U lab T , q lab ) = 4(m 2 Q + (q lab T ) 2 + (q lab L ) 2 - (q lab T ) 2 (U lab T ) 2 Cos 2 (φ lab ) (γ lab-cm T ) 2 ) (4.10)
where φ lab is the angle between the vectors U lab T and q lab and γ lab-cm T is the Lorentz relativistic factor calculated from just the transverse component of the 4-velocity

γ lab-cm T = 1 + (U lab T ) 2 .
From which one can derive

P lab T dP lab T = (s + U lab T ∂s ∂U T )U lab T dU lab T = (s -4( U lab T Cos(φ lab )q T (γ lab-cm T ) 2 ) 2 )U lab T dU lab T (4.11)
then defining the auxiliary quantity

J = (s + U T ∂s ∂U T ) = s -4( U lab T Cos(φ lab )q T (γ lab-cm T ) 2 ) 2 (4.12)
which allows to say that P T dP T = JUdU T = Jd 2 U T , then the normalization can be written as

d 2 U lab T dyd 3 q lab d 3 r lab JE lab Q Q × d 9 N Q Q d 2 P lab T dP lab L dq lab d 3 r lab = 1 (4.13)
In equation (4.13) we have expressed the integral of the normalization condition as a function of the variables in which we are interested in for studying the quarkonium production (y, U T ) and for applying our coalescence model ( q, r). Next, we express the 9 dimensional differential object as a function of quantities that are appropriate for the application of our coalescence principle.

d 9 N Q Q dyd 2 U lab T d 3 q lab d 3 r lab = JE Q Q d 9 N Q Q d 2 P lab T dP lab L d 3 q lab d 3 r lab = JE Q Q d 9 N Q Q d 3 P lab d 3 q lab d 3 r lab = JE Q Q d 9 N Q Q d 3 p lab Q d 3 p lab Q d 3 r lab (4.14)
Substituting equations (4.2) and (4.3) in equation (4.14) we obtain

d 9 N Q Q dyd 2 U lab T d 3 q lab d 3 r lab = JE Q Q × d 3 N Q Q d 3 r lab × d 3 N Q d 3 p lab Q × d 3 N Q d 3 p lab Q (4.15)
where the single quark momentum p lab Q, Q can be expressed as a function of the rapidity

d 3 N Q d 3 p lab Q = d 2 N Q d 2 p lab T,Q dN Q dp lab L,Q = 1 E lab Q d 2 N Q d 2 p lab T,Q dN Q dy Q . (4.16)
and the equation (4.16) For Q we find a similar expression. It is important to mention that the distributions

d 3 N Q Q d 3 r lab , d 2 NQ d 2 p T,Q , d 2 N Q d 2 p T, Q , d NQ dy Q and d N Q
dy Q are assumed to be normalized to unity, and are referring to the distribution of quarks (Q) and anti-quarks ( Q) at production.

Having obtained equation (4.15) and employing the result obtained in (A.57), we can derive then the probability density of the formation of Q Q = Φ bound state on according the by

dN i dy Φ = (2πh) 3 d 2 U lab T,Φ d 3 q lab d 3 r lab JE lab Q Q E lab Q E lab Q × d 2 NQ d 2 p T,Q d 2 N Q d 2 p T, Q d NQ dy Q d N Q dy Q × d 3 N Q Q d 3 r lab W Φ,NR (q cm , r cm ) (4.17)
where y Φ is the rapidity of the bound state Φ. Then knowing that d 3 q cm d 3 r cm = d 3 q lab d 3 r lab and

γ lab-cm E lab Q Q E lab Q E lab Q ≈ E cm Q Q E cm Q E cm Q .
This approximation is valid as long as the rapidity of the single quark Q and of the anti-quark Q are small. This is generally the case because our Wigner function has a peak for small relative momentum. We can express then

dN i dy Φ = (2πh) 3 d 2 U lab T,Φ d 3 q cm d 3 r cm JE cm Q Q γ lab-cm E cm Q E cm Q × d 2 NQ d 2 p T,Q d 2 N Q d 2 p T, Q d NQ dy Q d N Q dy Q × d 3 N Q Q d 3 r lab W Φ,NR (q cm , r cm ). (4.18)
On the other hand, for the distribution of the relative distance

d 3 N Q Q d 3
r lab between Q and Q, we can assume an expression of the form

d 3 N Q Q d 3 r lab = δ(r lab L )δ 2 (r lab T -r 0 ) (4.19)
where r lab L is the longitudinal component of the relative distance of the pair, while r lab T is the transverse component of the relative distance vector and finally r 0 is the initial relative distance between Q and Q when the Q Q pair is produced. This assumption is possible because in pp collision the formation of a bound state is a fast process which takes place immediately after the production of the Q Q pair. So, the distribution of the relative position will be a function of r 0 . Furthermore, as mentioned above, after formation the bound state will travel freely (no interaction) in the vacuum.

In longitudinal direction the component in the lab is related to that in the center of mass by r lab L = r cm L γ lab-cm , from which we can get δ(r lab L ) = γ lab-cm δ(r cm L ). After substituting the relation between the longitudinal components in the coalescence formula we obtain We can orient the initial separation r 0 along the Ox axis without losing generality, so we can write Then, knowing that r lab T,|| = r cm T,|| γ lab-cm T and r lab T,⊥ = r cm T,⊥ , where γ lab-cm T = 1 + (U lab T ) 2 , as well as d 3 r cm = dr cm L dy cm r dx cm r = dr cm L dr cm T,|| dr cm T,⊥ , we are ready to evaluate the integral over d 3 r cm in the coalescence equation (4.21), which takes the form

dN i dy Φ = (2πh) 3 d 2 U lab T,Φ d 3 q cm d 3 r cm JE cm Q Q E cm Q E cm Q × d 2 NQ d 2 p T,Q d 2 N Q d 2 p T, Q d NQ dy Q d N Q dy Q (4.20) × δ(r cm L )δ 2 (r lab T -r 0 )W Φ,NR (q cm ,
δ 2 (r lab T -r 0 ) = δ(y r )δ(x r -r 0 ). ( 4
d 3 r cm δ(r cm L )δ(y cm r )δ(x cm r -r 0 ) = dr cm L dr cm T,|| dr cm T,⊥ × δ( r cm T,|| γ lab-cm T Sin(θ) + r cm T,⊥ Cos(θ))δ( r cm T,|| γ lab-cm T Cos(θ) -r cm T,⊥ Sin(θ) -r 0 ) (4.25)
From the first two delta function we obtain the relations r cm L = 0 and r cm T,⊥ = -Tan(θ) γ lab-cm T r cm T,|| , then after substituting these relations in the third delta function we from which we derive that r cm T,|| = r 0 γ lab-cm T Cos(θ), leading us to r cm T,⊥ = -r 0 Sin(θ). Finally as we already relate the two components of the r cm T with r 0 , we can proceed now to see the dependence of the transverse relative distance vector with the initial separation

||r cm T || 2 = r cm T,|| 2 + r cm T,⊥ 2 = r 2 0 ((γ lab-cm T ) 2 cos 2 (θ) + sin 2 (θ)) = (1 + (U lab T ) 2 cos(θ))r 2 0
(4.27) As we already see that r cm L = 0, this mean that the Wigner function will be evaluated according with the equation (4.27). Thus, evaluating the r cm integral we have now

dN i dy Φ = (2πh) 3 d 2 U lab T,Φ d 3 q cm JE cm Q Qγ lab-cm T E cm Q E cm Q d 2 NQ d 2 p T,Q d 2 N Q d 2 p T, Q d NQ dy Q d N Q dy Q × W Φ,NR (||q cm || 2 , (1 + (U lab T ) 2 Cos(θ))r 2 0 ) (4.28)
Paying attention to the equation (4.28), we see that the angle θ appears only in the Wigner function dependence, and as we want to consider all possible orientations between the velocity U lab T and the pair relative distance direction, we will need to integrate over this angle. For doing so, we take the non relativistic expression for the Wigner (A.4)

W i,NR (q cm , r cm ) = 8g (πh) 3 e -r 2
dφW i,NR (q cm , r cm ) = 8g (πh) 3 e -(q cm σ) 2 h2 dφe -(1+u 2 T Cos(φ)) r 2 0 σ 2 (4.30) = 16g π 2 h3 e - (q cm σ) 2 h2 e -(1+ u 2 T 2 ) r 2 0 σ 2 I 0 ( u 2 T r 2 0 2σ 2 ) (4.31)
being I 0 especial Bessel function. Therefore, substituting the θ integral

dN i dy Φ = 16πg U T γ T dU T e -(1+ u 2 T 2 ) r 2 0 σ 2 I 0 ( U 2 T r 2 0 2σ 2 ) d 3 q cm JE cm Q Q E cm Q E cm Q × d 2 NQ d 2 p T,Q d 2 N Q d 2 p T, Q d NQ dy lab Q d N Q dy lab Q × e -( q cm σ h ) 2 (4.32)
Finally, the last step consist in expressing p Q and p Q in the lab frame as a function of q cm r and U T , and then substituting those relation in the single quark production densities

d 2 NQ d 2 p T,Q d 2 N Q d 2 p T, Q d NQ dy lab Q d N Q dy lab Q .
Implementing the component of the q T parallel and perpendicular to the vector U T

q lab T,|| = γ lab-cm T q cm T,||
we can write the expression for the single quark(anti-quark) momentum in the lab frame as following

p lab Q = ( √ s 2
U T + γ lab-cm T q cm T,|| ) ÛT + q cm T,⊥ ÛT,⊥ + q cm z e z (4.33)

p lab Q = ( √ s 2
U Tγ lab-cm T q cm T,⊥ ) ÛTq cm T,⊥ ÛT,⊥q cm z e z (4.34)

where ÛT,⊥ is a unitary vector in a transverse direction with respect to U T , while ÛT is an unitary vector in the direction of U T , √ s 2 = m 2 Q + q cm,2 and γ lab-cm T q T,|| is the component of the relative transverse momentum q T which is parallel to U T and thus get affected by the γ lab-cm T . It is also possible to obtain an expression for factor J as a function of the components of the vector q T from the equation (4.12). First, it is important to note that

Cos(θ lab )q lab T = q lab T,|| = γ lab-cm T q cm T,||
Knowing that, one gets after substituting the components of the vector q cm T in (4.12)

J = s -4( U T γ lab-cm T q cm T,|| (γ lab-cm T ) 2 ) 2 = s -4 (U T q cm T,|| ) 2 γ 2 T = 4(m 2 Q + q 2 z + q 2 T,⊥ + q 2 T,|| γ 2 T ) (4.35)
Similarly, the single quark (anti-quark) transverse momentum can be expressed

Chapter 4. Charmonium Production in pp Collisions

as a function of the transverse momentum of the center of mass (P T ), the relative transverse momentum q T and the angle θ as follows p T,Q = (P T ) 2 + (q T ) 2 + P T q T Cos(θ) (4.36)

p T, Q = (P T ) 2 + (q T ) 2 -P T q T Cos(θ) (4.37)
The final step in our analysis is integral over the angle θ, the z-component of the relative momentum q cm z , and over q cm T,⊥ in that specific order. Integration over U T can also be carried out as the final step of this strategy if it is desired to study quarkonium production only as a function of rapidity ( dN dy ). But the study of quarkonium production with respect to U T (P T ) is also of great interest for the present work.

That said, if in equation (4.32) we assume that there is no correlation between the rapidity of the quarks (anti-quark) and the position and momentum variables referring to the heavy quarks pair, the equation (4.32) can be rearranged in the form

dN i dy Φ ≈ P Φ × dN(y Φ ) dy Q dN(y Φ ) dy Q (4.38)
where in the equation (4.38) a factorization between the single quark (anti-quark) rapidity and the pair related coordinates is assumed in base of the uncorrelated assumption between this two group (single quark rapidity and pair relative coordinate) of variables as well as the assumption that

P Φ = 16πg U T γ T dU T e -(1+ U 2 T 2 ) r 2 0 σ 2 I 0 ( U 2 T r 2 0 2σ 2 ) d 3 q cm JE cm Q Q E cm Q E cm Q d 2 NQ d 2 p T,Q d 2 N Q d 2 p T, Q e -( q cm σ h ) 2
(4.39) the later equation, provide us with the probability that a single Q Q end up forming a quarkonium state Φ, which from now on will be recalled as conditional probability. On the other hand, the equation (4.38) give us the total probability of quarkonium state formation (Φ) by factorizing the single pair probability and the number of heavy quark produce in a pp collision event.

Taking for instance, the production of the J/ψ in pp collision at √ s = From this fit, one derives the distribution of q T , by substituting the relations of the equation (4.37) and then integrating the distribution product d 2 Nc

d 2 p T,c d 2 Nc
d 2 p T, c on θ and p T . This procedure lead us to

q 2 T = 4.35GeV/c (4.41)
Hence, one has

q 2 T × σ r ≥ n D × h 2 (4.42)
where n D is the number of dimension of the integration differential (n D = 2) and from where we can derive that σ r ≥ σ r,min = h 4.35 ≈ 0.044 fm the average transverse relative distance (spatial width). As for r 0 (which represent the initial separation of the pair members), one can take two approach to tackle down the lower bound limit issue of r 0 , the first A: is simply identify the initial separation to just equal to the average relative transverse separation ( ||r 0 || 2 = σ r ) of the Heisenberg principle. Meanwhile the approach B: introduce a distribution of ||r 0 || with some r.m.s σ r for example

P(||r 0 ||) = 2||r 0 || σ 2 r e - ||r 0 || 2 σ 2 r (4.43)
in a similar way that was discussed in [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF], in which the value for σ r was estimated as function of mass of charm quark

σ r = 8 3 h m c ≈ 0.2 (4.44)
The analysis of the approach A is rather trivial, in just enough to determinate (or chose) a proper value for the average transverse separation σ r and therefor of r 0 , meanwhile the approach B, is a bit more complicated. In order to applied the approach B, one needs consider also the distribution (4.43) at the moment of computing P φ . Then, including the distribution (4.43) inside the conditional

P Φ = 16πg U T γ T dU T I B (σ r , σ, U T ) d 3 q cm JE cm QQ E cm Q E cm Q d 2 NQ d 2 p T,Q d 2 NQ d 2 p T,Q e -( q cm σ h ) 2 (4.45)
with

I B (σ r , σ, U T ) = dr 0 e -(1+ U 2 T 2 )
r 2 0

σ 2 I 0 ( U 2 T r 2 0 2σ 2 ) × P(r 0 ) (4.46) = 2 σ r ∞ 0 dr 0 r 0 e -( 1+ U T 2 σ 2 + 1 σ 2 r ) I 0 ( U 2 T r 2 0 2σ 2 ) (4.47)
being the integral over all possible values for the initial separation r 0 , one has for this particular integral

∞ 0 dxxe -βx 2 I 0 (αx 2 ) = 1 2 β 2 -α 2 (4.48)
from where we get substituting the previous integral

I B (σ r , σ, ...) = 1 (1 + σ 2 r σ 2 ) × (1 + (1 + U 2 T ) σ 2 r σ 2 ) (4.49)
Taking m c = 1.5 GeV (charm quark), one can obtain the dP J/Ψ du T "spectra" within approach A (||r 0 || = σ r ). Employing two different values for σ (σ = 0.303 fm) corresponding to the value adopting in [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF], while σ = 0.35 fm is an optimize value obtained from the internal notes summary analysis. And for the choice made for the initial separation r 0 values, the selection seems to be more difficult. The first interesting value to check is r 0 = 0 fm, although this is a clear violation of the quantum mechanics laws, secondly r 0 = 0.044 fm, is the smallest value allowed by Heisenberg principle applied over relative momentum distribution, and finally r 0 = 0.3 fm that was the standard prescription in EPOSHQ. , especially when r 0 ≥ σ r = 0.044 fm. This can be explained from the point of view of the range of the Wigner distribution peak, if the quark-antiquark pair initial separation is very large, they have less chance to be "caught " by the distribution and therefore less probability to form a quarkonium state. On the other hand, the choice of σ ≥ r 0 the influences over conditional probability spectra are also very notable. That is, for larger values of σ implies smaller width in momentum space and therefore less chance for coalescence to occur together with the fact that right after creating the HQ pairs will be immediately within the range of the Wigner function peak. In the following figure, the same differential quantity is presented but now employing the approach B (a distribution according equation (4.43) for the initial separation r 0 . Even for this approach, the same trends are observed, showing us that the features before mentioned are quite general and does not change even if we change the production separation (r 0 ) within range of values allowed by the quantum mechanics. As we can see, while increasing the value of σ r the conditional probability decrease, this is due to the fact that larger values of the distribution width, allows larger value of r 0 , which could leads us to cases in which the HQ pair is outside the Wigner function peak. Now, comparing both approach (A and B) chosen the values of r 0 in the case of approach A to be equal to r.m.s of the distribution (4.43) σ r value in the approach B, we can notice that for small values of r 0 (σ) the approach B offers higher probability production. This is due to the fact in the B approach, as we are using a distribution for r 0 with r.m.s equal to σ r , there still always some chance of obtain small values r 0 which leads to some higher probabilities cases.The later, is more notable for larger values of r 0 (σ r ), because for those cases, the difference between the selected value of r 0 for the approach A and the value of r 0 coming from the distribution (4.43) in the approach B can be more large.

To complete the conditional probability analysis, the values for P J/ψ (integral of the dP J/ψ dU T

) for the approach A and B are show in the following two tables. From the general considerations from the analysis of the results of the conditional probability, one might conclude that the later (conditional probability) provide important insights quarkonium formation. We can understand (by consider only a single Q Q) the impact of some quantities such like initial separation r 0 over the probability for quarkonium production. As the results show the influence of this quantity is important and for that the selection of the proper values for those "free" parameters is critical, but at the same time our model allow us to evaluate quantitatively the effect over the probably and thus initiate some discussion about it. Although, there are several good estimation for the values of those parameters [START_REF]Heavy-quark production in heavy-ion collisions[END_REF], there are not too much further discussion on this regard.

p T spectrum of quarkonium production in pp collisions

Moving forward in our analysis, we are now ready to to evaluate the inclusive J/ψ production in pp events. The relationship between the inclusive differential J/ψ production and the conditional probability can be expressed The black dotted line in the figures indicates the experimental reference obtained by considered the experimental value for dσ c dy ≈ 1 mb [119] reported at mid-rapidity and the experimental value for dσ J/ψ dy = 6µb reported in [118] from which we derive a ratio 0.006 which lies just in the bulk of our predictions. It is important mention that we neglect the contributions from the feed-down (that should be take into consideration in order to make a rigorous comparison with experimental inclusive production data).

Focusing, from now on, only in the approach B, we present now a comparison of the P T distribution prediction of our model with experimental data published in [118]. But before that, we need to obtain the explicit expression for the quantity that have been already plotted in the previous section. Then, comparing the results obtained from the application of our model (equation (4.56) to the mid-rapidity J/ψ production data at 5.02 TeV in pp collision provide by the ALICE collaboration. In the following figure the comparison between our prediction and the experimental data is presented It is important to recall that the calculations from equation (4.56) do not include the feed-down from the excited states of charmonium or B-mesons. The later (B-mesons feed-down) being significantly large especially at higher values of p T , which can explain why our model underestimate the J/ψ production at large p T . Even so, despite the lack of the feed-down contributions our model reproduce rather well the experimental data shape all the way up to p T = 5 GeV, being (σ =0.35 fm and σ r = 0.2 fm) the optimal values for the "free" parameters to do so. These two values will become the reference values for future calculations in the present work, especially in the next chapter.

Chapter 5

Quarkonium production in AA collisions

As explained in chapter 3, Remler model consists fundamentally in considering the contributions of all possible pairs between the particles that make up a composite state (mesons or light nuclear fragment) and calculating the probability of formation of those states for all possible pairs in a given time. This is done through the projection of the phase space density in relative coordinates of the final state (composite state) onto the initial state (the separate components). As discussed also in chapter 3 a J/ψ can only be formed below a given dissociation temperature (T diss ) but it can be formed from all c and c independent of the vertices from where they originate. The later is the first major difference between AA collisions and pp collisions. In pp collisions we rarely have creation of a J/ψ from one c c pair consisting in a c-quark from one vertex and c antiquark from a different vertex. Another major difference is that in AA collisions the formation of bound state is delayed with respect to pp collision case due to the medium temperature, meaning that c-(anti)quarks need to arrive or wait for the local medium temperature to cool enough to permit the bound state formation.

From now on, we will name "diagonal contributions" those J/ψ's made from heavy quarks from the same vertex, making an analogy with a matrix arrangement in which we can associate the indices of the columns with the number of quarks Q i and the rows with the index of the anti-quarks Qj . In this way, the diagonal elements will be those that come from the same "mother parton" and therefore will have the same associated index. In AA collisions the production vertex of c c pair are located at a relatively short distance. This conditions the existence of new cases in which a c from vertex 1 can be relatively close to a c from vertex 2, and therefore a reasonably large probability of creating a quarkonium state can be expected. These contributions will be intuitively called "off-diagonal contributions".

Another big difference between the physical scenarios of pp and AA collisions is the influence of the medium (QGP) on the formation of the quarkonia. The discussions about how the influence of the environment was taken into account in our model is quite complex and will be divided into several parts for a better understanding. First, in the equation (3.36), the initial (primordial) probability P prim depends on time due to the dissociation temperature for each quarkonium state T i diss and the position of the heavy quarks inside the plasma (local temperature). This is mainly since not all c-(anti)quarks arrive at this threshold at the same time and are only "available" at later times for the formation of bound states. The formation of bound states in a dense and hot medium will be modified with respect to that in the vacuum due to the screening effect on the cc interaction potential which modifies the Wigner density function of the bound state. The latter comes from two effects:

1. emergence of an additional rate term (local rate) as described in section (3.4);

modification of two-body density function (Wigner function) through the

implementation of a width (σ = σ(T)) which depends on the temperature of the medium.

An essential ingredient of the present contribution is to have a well-calibrated energy loss model. As is the consensus in the HIC community, the medium induced radiation energy loss increases their relevance with increasing energy (p T ) of the c-(anti)quarks. This is very important in the presented work because this type of energy loss was neglected in our the calculations that will be discussed here. However, elastic energy loss was correctly calibrated as the D meson R AA is correctly reproduced [START_REF] Aichelin | Heavy-flavor azimuthal correlations of D mesons[END_REF][START_REF] Nahrgang | Heavy-flavor azimuthal correlations of D mesons[END_REF]. Therefore, one can have some degree of confidence modelling the heavy quarks interactions with the medium only with collisional energy loss.

All the aforementioned considerations are the main ingredients of our model, so in this chapter, I will present the results referring to the direct J/ψ production, and I will discuss the main features and consequences of our model and the aforementioned considerations.

Proton-Proton scaled production

One of the main objectives of our models is to describe the suppression of the quarkonium (charmonium) production in AA collisions. The first step towards this goal is to study what our model predicts for the J/ψ production in AA collisions immediately after the c c pairs have been produced in elementary binary collisions. This calculation will be performed before the QGP pre-equilibrium phase ends (τ ∼ 0.01 fm/c) and therefore before the medium can affect the production. This implies that in this calculation we assumed that our c-quarks are in the vacuum (σ vac = 0.35 fm).

If one wants to estimate pp-scaled production (pp production already multiplied by the number binary nucleon-nucleon collision) directly from an AA collision, we need to limit our study to those pairs created from c and c that come from the same production vertex (diagonal contribution). Therefore, the considerable enhancement observed in [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF] coming from the contributions of pairs formed by c and c coming from different production vertices (off-diagonal contribution), will not be present in this result.

The pp-scaled production probability of charmonium in our model is given as

P Φ i pp = N diag ∑ k=1 W NR,Φ i (r cm rel (k), p cm rel (k)) (5.1)
where N diag stands for the number of pairs coming from the same production vertex (one vertex per pair), while r cm rel (k) ans p cm rel (k) are the relative distance and relative momentum associated to the k th pair. This production will be our reference for the calculation of nuclear modification factor R AA . In figure 5.1, after comparing the results of the application of our coalescence model to genuine proton-proton collisions and the results obtained through the selection of diagonal pairs in PbPb collisions, the results are extremely similar; especially since CNM effects such as (anti-)shadowing have not been considered.

With the relativistic corrections applied to the original Remler formalism, our model can deal in a very satisfactory way with the relativistic quarkonium (J/ψ) formation as long as the relative momentum between c and c that form the bound state is small (p r < m Q ) and for p J/ψ T 5 GeV/c. On the other hand, our model underestimates the charmonium production at large p T reported by the ALICE collaboration [119]. This feature could be due to the fact that, in the present stage our model, the underlying c c distribution on which our coalescence mechanism takes place corresponds to uncorrelated c and c 1-body distributions and does not contain gluon-splitting contributions which are known to play an important role at large p J/ψ T It is also important to emphasize (as was commented in chapter 4) that our model calculates the direct J/ψ production. Therefore, one has to be cautious to compare 5.1. Proton-Proton scaled production 115 our results with the experimental prompt data since we only compute the direct J/ψ production component.

Another interesting perspective to evaluate in our model is to study the production of J/ψ in the forward rapidity range. The study of this case will be restricted in the present work to the case of pp-scaled production only. The reason for this decision is the way in which the time-evolution of c-quarks in the plasma is treated. As discussed in Chapter 3, the temporal evolution of c-quarks is carried out through a constant Minkowski time-step. This implementation of the dynamic evolution of the quarkonium formation probability cannot be trivially extrapolated for the case of forward rapidity (2.5 < y < 4.0). First, the relatively large value of the rapidity would imply that the time step value needs to be increased properly. Here, the time at which the calculation was performed was also τ ∼ 0.01 f m/c, as in the case of mid-rapidity. This is a difficulty due to the time dilation in the case of forward rapidity with respect to the case of mid-rapidity. This means that the timescale used may not be suitable for forward rapidity (in particular, too early).

Still, our predictions for pp-scaled production at forward rapidity (5.2) are good , showing that future applications of our model for forward rapidity are possible. Even so, the best way to perform this calculation is through one of the options offered by EPOSHQ, which allows us to move our reference system from the center of mass of the collision to a new reference system that fulfills that y=3. [START_REF] Albacete | Gluon saturation and initial conditions for relativistic heavy ion collisions[END_REF]. By doing this (placing ourselves at the center point of the forward rapidity interval), we can repeat our analysis implemented for the case of mid-rapidity. This in turn, will allow taking into account the time dilation effects inherent to the particles in forward rapidity values.

Primordial contribution

The primordial production refers to the formation in the early phase of the evolution of QGP (or even before the formation of plasma) of bound states (charmonium) due to the closeness between the c c pair numbers produced in the initial hard scattering.

As discussed in chapter 3, in our model, the contributions of each c ( c) to the primordial production are calculated once they arrive at the production regions (regions in the plasma whose local temperature is T < T diss ). The latter, together with the use of a Wigner function that depends on the local temperature (local basis), represents a difference with the predictions reported in [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF], which were obtained by using a Wigner function that does not depend on the temperature (vacuum basis). Hence a direct comparison between our results and the predictions in [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF], in which an initial enhancement in AA production with respect to the equivalent pp was found, is not possible for the results presented in this work. However, in the initial stages of this thesis, the vacuum basis was used to establish the preliminary predictions.

Definitions

The formation and thermalization of the QGP is one of the most debated topics in the community. The consensus is that after a small time window after the initiation of the nucleus-nucleus collision [0.1-1] fm/c [START_REF]Heavy-Ion Collision-hot QCD in lab[END_REF][START_REF] Eskola | Formation and Evolution of Quark-Gluon Plasma at RHIC and LHC[END_REF] the formation and thermalization of the QGP has been completed. As already discussed in this work, c c pairs are produced in the earliest stages of plasma and then travel in the medium carrying out many interactions with plasma particles. Due to the high temperature in the early phase of QGP, the formation of bound states (quarkonium) is unlikely because the temperature of the medium is too high (T > T i diss ), see also (table 2.1 in section 2.2.2).

These temperatures above T J/ψ diss = 0.4 GeV (for the case of J/ψ) prevent the stability of any bound state due to the high rate of dissociative processes such as

Φ + g → Q + Q
and due to the strong screening of the binding potential. In our Monte Carlo simulations, the primordial production is thus a dynamic process, that takes place over a relatively long time span. In other words, there is no sudden production at a given time point. Each individual c c pair will be considered for its contribution to the primordial production once both these quarks have evolved long enough in the QGP to be located in a region whose temperature is lower than the dissociation temperature of the quarkonium state T i diss while still not being part of any final bound state like D or D mesons.

The c-quarks are produced in different elementary collisions distributed throughout the overlap region of the nucleus-nucleus collision. Therefore the different c-quarks will fulfill at different times the condition to be "available" to form a bound state, due to the anisotropy in the energy density and temperature in the QGP. Concretely, the following process is adopted: When a c-quark arrives for the first time in a region with T < T diss we calculate its probability to form a J/ψ with all other c which are already satisfying this condition. The sum of all these contributions for all c-quarks is what we call primordial production (or contribution). On a calculation level, one introduces the quantity

P Φ i ,prim k (t c ) dy k = n c ∑ l=1 δ(y k -y Φ )W NR,Φ i (r cm (k, l), q cm (k, l)) (5.2)
pertaining to the k th c-quark, where n c is the number of anti-quarks that are "active" (that are inside the QGP, in a region of the plasma which a local temperature lower than T i diss threshold and also are single quark) at the instant of time t. The number of Chapter 5. Quarkonium production in AA collisions primordial J/ψ created between time interval [t, dt] is then defined as the sum

dN prim dy (t, t + dt) = N first (t,t+dt) ∑ k=1 P Φ i ,prim k (t c )) dy k (5.3)
where N first (t, t + dt) stands for the number of c-quarks which fulfill the temperature condition (T(x c ) < T i diss ) for a first time t c taking place in this time interval [t, t + dt].

Results for the primordial production

The first result of interest turns out to be the temporal evolution of the primordial production described by the equation (5.3). It gives the total probability of charmonium formation "at" a given instant of time (t) An exciting aspect of these results is that they reflect the average time scale that the c-quarks are ready for forming bound states.

As shown in the figure 5.3, most of the c-quarks arrive in regions of the plasma in which the formation of bound states (those regions whose temperature is lower than the dissociation temperature T diss ) is possible relatively early (0.5 ∼ 2.5 fm/c). This result is quite expected because, as we know, c-quarks are produced in the early phase of a HIC. Therefore, the sporadic cases in which they arrive in areas where quarkonium (charmonium) production is possible relatively late (t ≥3 fm/c) correspond to quarks created in regions of the very hot medium, where they need more time to cool down or due to those c-(anti)quarks created at the boundaries of the rapidity range selected for this figure (|y| < 0.9) which undergo some time dilatation.

One can also notice in inter sting effect when the time evolution of production in the different centrality bins is compared. The primordial production for the case of very central collisions (for example, 0-20%), in addition to the higher production due to the larger initial number of c c pairs, is also more retarded compared to the less central cases. This is because the more central the collision is, the hotter and denser the plasma produced is, which causes that, on average, the c and c reach the production regions later. Another aspect shown in the results of figure 5.3 is the similarity between the cases with interaction potential ON and OFF. This last aspect can be studied in more detail while looking to the p T -spectra of the primordial time-integrated (see figure 5.4) production. The best way to get a complete picture of the influence of the medium (elastic collision), which is the only energy loss mechanism considered in the evolution of the c( c) and of the effects of the potential, is to study the behavior of the production probability when we activate or deactivate elastic collisions with the medium on the one side and potential interaction between c and c on the other side. From that analysis, four possible scenarios were studied for the evolution of c and c, based on considering the combinations of c c potential interaction (ON and OFF) and medium elastic collision (ON and OFF).

From the analysis of the results in the figure 5.4 and focusing first on the effects of elastic collisions of the c-(anti)quarks (Med Coll) with medium particles on the primordial production, we see that when elastic collisions are active, the p T -spectra is shifted towards lower p T : Those c and c produced initially with higher p T (energy) lose energy through collisions with bulk particles. Although the primordial production happens mostly early (see figure 5.3) in the QGP evolution, the heavy quarks have already suffered enough collisions to be sensitive to this shifting effect. This indicates that the frequency of elastic collisions of the c-(anti)quarks at the early stage is very high. On the other hand, the effects of the c c potential interaction at the stage of primordial production are almost negligible. This can be explained by the strong screening of the potential. As stated earlier, at the primordial stage, the temperature is very high; the screening effects nullify the action of the potential; therefore, the primordial production is mainly dictated by the fact that the c and c quarks are still relatively close at that stage. To conclude this section on primordial production and offer a necessary reference level on the general production of c-(anti)quarks in the different centrality bins. In the following table (5.1) the characteristic values of the number of nucleon-nucleon binary collisions N coll , the average number c c pair ( dN c c dy y=0

) and initially produced per unity of rapidity, and finally the total average number (

dN prim J/ψ dy y=0
) per unity of rapidity and the ratio between the primordial J/ψ and the number of produced c c per unity of rapidity (

dN prim J/ψ dy / dN c c dy y=0
) .

Analyzing the results presented in the table (5.1), we can see that the number of c c pairs per unit of rapidity scales as expected with N coll . On the other side, the production of primordial J/ψ does not scale with respect to the number of c c pairs produced initially in the same way as that exposed in [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF]. According to what is stated in [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF], one should expect an almost quadratic scaling of the number of primordial J/ψ with respect to the number of c c pairs per unity of rapidity, meaning that the ratio

dN prim J/ψ dy / dN cc dy y=0
should increase linearly as we move from peripheral to a central collisions. In our results however, we obtain a different trend (see 5 th column of table 5.1). The essential difference between our results and the predictions made in [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF] lies in the way in which the primordial production is calculated. In [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF] the primordial production is calculated at a very early stage of the collision, right away after all c c pairs are produced and employing a Wigner function which corresponds to a vacuum basis and thus to does not take into account the dissociation temperature. However in our case, the primordial production was obtained through a Wigner function that depends on the local temperature of the medium (local basis) and is distributed in time (delayed production). To see in more detail the influence of centrality beyond the differences derived from the initial number of c c pairs produced, the temporal evolution of the primordial production normalized by the square number of c c per unity of rapidity ( dN cc dy ) 2 for the case of interaction potential and medium elastic collision (ON) is presented in figure 5.5. In this figure, we can see that for very early times [0-1.0 fm/c], the primordial production per c c pair is more suppressed as a consequence of the hotter and denser QGP created in this very central collisions. However, as we move forward in time and more c and c arrive in cooler regions, the most central collision starts to dominate the production due to the large number per c c pairs, which induces non-linear effects in the production probability per c c pairs coming especially for the numerous off-diagonal contribution (the same effect as the one investigated in [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF]).

Analysis of heavy quarks trajectories

One of the essential ingredients in our model is the semi-classical trajectories of the HQ, which reflects the influence of the medium on the quarks as they propagate through the plasma and the interaction potential that we are considering.

An interesting way to study the consequences of the four scenarios presented in

Analysis of heavy quarks trajectories

123

the previous section (interaction potential: ON and OFF as well as medium elastic collision: ON and OFF) on the trajectories of the c (and c) quarks is to pay attention to the behavior of the relative coordinates (r rel and p rel ) as a function of p T . In particular, one can study the total number (integrated over time) of c c pairs that are closer than a specific "distance" in phase space, at least at a given time instant (see figure 5.6). From figure 3.3, we know that the interaction potential keeps the c and c quarks close to each other for longer times, which implies that the average value of the relative distance r cm rel is lower than the case when we do not consider the potential. On the other hand, from our previous results (figure 5.4), we know that the elastic collisions with the medium shift the population of HQ towards low p T , which implies that the HQs lose energy (momentum). As it has less energy (momentum), the average relative momentum tends to be lower than the cases where medium elastic collisions are not considered.

According to this analysis, one should expect that the scenario in which both (interaction potential and medium elastic collision) are ON will be the most favored in terms of production of close c c pairs which will be, besides, located at smaller p T . This trend is clearly observed on figure 5.6 (green curve). Product to combine two effects that help in the production (see equation A.4 ) of bound states: 1-the decrease of the kinetic energy of the quarks product of the medium elastic collision, which helps to decrease the average relative momentum 2-the interaction potential which helps c c pairs stay relatively close for longer, thus decreasing the average value of the relative distance of the pairs. This also justifies why the case in which both (Interaction potential and medium elastic collision) are OFF is less favored (yellow curve on figure 5.6).

On the other hand, in the "crossed" scenarios (Interaction potential (ON) and medium elastic collision (OFF) and vice versa), the behavior with respect to p T is defined by the dominance of one of the effects mentioned above while the other is absent. When there is only medium elastic collision (no interaction potential), the shifting effects makes this scenario low p T stronger than the scenario in which there is only interaction potential (no medium elastic collision). While at high p T we see an inversion of the situation due to the absence of medium elastic collision in the only potential case (blue curve of figure 5.6)

Local Rate

As commented in chapter 3 (section 3.4), the local rate is a natural consequence (A.62) while readdressing the Remler model results and considering a temperature (time) dependent Wigner density function. The influence of the medium, more specifically the evolution of the local temperature in the region in which the c c pairs are located, alters the wavefunction (and therefore the density operator) that describes each of the quarkonium states (charmonium), due to the modification of the charmonium's Hamiltonian.

During the implementation of the local rate computation, the numerical precision is limited by the discretization of the time derivative of Wigner width (σ(t)) (A.63) and the definition of the medium local temperature where the c c pair center of mass is located. This is because if the relative distance is large, the local temperatures associated with region in which the c (T c ) is located and the one in which c (T c) is located are different 1 . The simplest and most direct way to assign a temperature to the two-body system (pair) is via the "average" temperature between the two regions where the c and c are located.

T c c = 1 2 (T c + T c) (5.4)

Local Rate 125

Once the temperature associated with the pair as a whole has been established, the temporal derivative of the Wigner function width (σ) will be evaluated in our simulations, employing the approximation

dσ(t) dt ≈ ∆σ(t) ∆t (5.5)
where ∆t is the computational time step measured in the pair's center of mass. To begin the study of the properties of the local rate, I first present the p T -spectra of the integral local rate.

Integrated Rate = ∞ t first Γ c c loc (t )dt (5.6)
where t first is the time instant in which the latter of the two (c or c) contributed to the primordial probability, that is first time from which the c c pair is forming a bound J/ψ. It is important to note that for the trivial case of c c with a finite center of mass velocity, the time derivative of the Wigner function ( time derivative of the Gaussian width σ) is taken with respect to the computational time measured in the lab (see equation 3.42). But when we calculate the integral over the evolution time, the integration variable dt is also taken from the computational frame. So both "dt" cancel each other out, including the gamma factors that come out when one calculates the value of the time integral local rate in the center of mass frame. By moving through different centrality bins (figure 5.7), we are altering the number of participants in the AA collision. This influences directly the initial parameters of the QGP (its density and temperature distributions as well as the geometry of the collision). Generally, as the centrality bins increase (larger centrality percent), the AA collision will produce a less dense and less hot plasma, which cools faster, and the initial production of particles (including heavy quarks) is considerably reduced. As was observed in the case of primordial production, the contribution of the integrated local rate to the total probability also decreases with increasing centrality.

Contrary to the case of primordial production (particularly in the case of integrated production over time), the influence of the interaction potential is notable for the local rate, as shown in the figure 5.8. The contributions to the local rate come from a more extended time range (from the primordial stage to the end of the QGP evolution), and for later times the c c pair interaction potential becomes stronger. So those c c that are close in phase space will undergo a strong modification in their dynamic evolution due to the interaction potential which becomes stronger as the medium temperature decrease. In particular, those trapped in the interaction potential well, will hardly escape unless a very hard gluon breaks the entanglement. This mean that the that generic effect of the potential is to increase the local rate and thus to produce more final J/ψ As for the role of the elastic collisions with the medium, in the very low p T sector, the elastic collisions with the plasma particles establish some hierarchy of production. This can be explained because the c-(anti)quarks in the low p T sector have low kinetic energy (especially in mid-rapidity), combined with the high population of plasma particles (light quarks and gluon) in the low p T sector. This leads to an extremely high rate of collisions, forcing the c-(anti)quarks in this sector to a Brownian-like dynamic. Therefore the heavy quarks pairs stay longer close to each other, boosting the probability of quarkonium formation. According to the definition of the local rate and of the Wigner function -eq. (A.63) and eq. (A.4) -that define the probability of formation of a bound state, the pairs with relatively small values of p cm rel and r cm rel are the most favored to form the bound states. The latter implies that the case in which the interaction potential and medium elastic collision are active, the integrated local rate gets favored from reducing the average relative momentum due to the medium elastic collision and decreasing the average relative distance due to the interaction potential (see also 5.6). In the same way, the opposite case, Med Coll and Potential (OFF) will be the most affected by these effects. In that respect the local rate hierarchy follows quite nicely the hierarchy alreay observed in figure 5.6. The other two scenarios ( the "crossed" ones) also inherit the behavior in tune with the results presented in that figure.

Collision rate

The collisional rate contribution to the total probability (in other words, the rate introduced by Remler and described in section 3.1 -3.3) is due to the change in the probability of quarkonium state formation by a collision with a third particle from the medium or, as we discussed later on, due to the interaction potential.

This contribution is strongly dependent on the frequency and the physical nature of the interactions that heavy quarks undergo during their movement in the plasma. We assume that the mean free path is sufficiently long that before and after the interaction (c c potential or collision), the c-(anti)quarks has well-defined phase space variables (r c, c, p c, c), which allow to compute the c c pair projection on a given charmonium state (J/ψ)

Collision criteria

In order to analyse the production of Quarkonia in our approach, in which the position and momenta of all particles are updated with a finite time step ∆t we introduce the collision coefficient. In the derivation of the collision rate, (A.4), it was established that that the probability to find a quarkonium changes only if one of the pair members suffered from a collision. However, employing a potential interaction as was described in section 3.5 also modifies the momentum during this time step ∆t (however mostly very weakly) what could be spuriously confused with a genuine collision by the Remler algorithm. In order to separate both contributions, we introduce the collision coefficient

η = |p c, c n -p c, c n-1 | |p c, c n | (5.7)
where p c, c n and p c, c n-1 are the momentum of the c-(anti)quarks (c and c) at the time step n and n -1 respectively. We consider that only c-(anti)quarks that undergo such a momentum variation larger than a given η cut contributes to the collision rate. The choice of the value for the collision coefficient is crucial because it can directly discriminate the effect of many processes that the c ( c)-quarks suffered during their evolution. On the one hand, choosing too large value can lead to the rejection of genuine collisions with the medium. On the other hand, choosing a value very close to zero can create unwanted effects like uncertainties due to interpolating the quarks momentum at a given time or any uncertainty derived from the time step of the evolution.

To begin the analysis of the effects due to the modifications in the c( c) momentum in the production probability, I will present the results corresponding to the p T -spectra of the average variation in the Wigner density function ( figure 5.9) due to a single modification, defined as

∆W J/ψ coll (p T ) = ∑ N c, c coll i=1 ∆W J/ψ i (p T ) N c, c coll (p T ) (5.8)
where N c, c coll (p T ) stands for the total number of collision suffer by the c or c,and finally ∆W J/ψ i is the variation in the probability of a J/ψ formation from a given pair if one of the members suffers a collision.

The modification were considered independently of whether they come from c c interaction potential or medium elastic collision. We see in the results of the figure 5.9, that for very low and high p T the variation of the probability is very small while around 4 GeV/c, this distribution finds its maximum, decreasing again as we move towards larger p T . This behavior can be explained analyzing the conditions under which we can obtain a large variation in probability. To produce relatively large variations in probability, the relative coordinates (mainly the relative momentum) must vary significantly between the initial and final instant, and at the same time the values of the relative coordinates should be relatively small for at least one of the two momenta (initial or final), which will cause the Wigner density value to be high.

Based on this analysis, we can conclude that in the low p T sector, while the values of the Wigner function are relatively high, the collisions in this sector are very soft and produce very small variations in relative momentum. So the difference between the final and initial values of the Wigner function is relatively small.

On the other hand, in the high p T sector, collisions with medium particles can produce relatively large variations in relative momentum, but the values of the Wigner function are relatively small compared to the values in the low p T sector. This will cause the absolute value of the variation in the Wigner function also to be small.

As complementary study is presented in figure 5.10. Here we display the integrated collision rate obtained for 0 -20% centrality for elastic collisions with the medium and c c potential interaction active for different collision coefficients ( η = 0, 5 and 10%). As long as we consider the collision coefficient equal to zero (η cut = 0%), any modification to the total momentum of the c ( c)-quarks from one time step to the next will be considered as a collision. Of course, as mentioned before, variations of the total momentum due to the interaction potential and those derived from the interpolation process will be taken into account for this case. Although the variations in momentum due to the interaction potential are a consequence of the physics of our problem and are welcome in our analysis, the variations derived from other sources (uncertainties) need to be removed therefore η cut = 0% is not an ideal value.

In figure 5.11, we present separately the momentum variation stemming from both processes considered separately. From this figure, one concludes that the modifications of the momentum of the c and c due to elastic collisions with medium particles (blue curve in the figure) are dominant. The sum of the effects of these collisions (they combine recombination and dissociation processes), are negative (figure 5.10) or in other words and according to the equation (A.4) the value of the probability of formation of a charmonium state after a collision is less than the probability of formation before. For the part of the elastic collisions with the medium's particles, the dissociation processes are dominant over the recombination processes. The latter is important when analyzing the behavior of the collisional rate for finite values of η cut . By increasing the value of η cut (5% and 10%) we neglect more and more the change of the probability due to potential interaction which took place in the time step which we have considered because in general the momentum change due to potential interactions is small and concentrated in the low percent part of the graphic 5.11. The interaction potential is always attractive and tend to increase the J/ψ production probability in a given time step. At the same time, medium elastic collisions always dominate the number of interactions (collisions), which implies that as we rise η cut we will also eliminate a large number of momentum variation due to medium elastic collision than those due to c c potential interaction. This yields to a decrease in the absolute value of the collision rate, because while discriminating more medium elastic collision processes and therefore removing from our analysis recombination and dissociation as well, the dissociation processes are more dominant, meaning that they will suffer more from this discrimination.

Influence of the Time Step

As our model aims to obtain a dynamic description (temporal evolution) of the formation of quarkonium states, the value of the time step is crucial. As already said, the collision rate is fed by two sources: elastic collisions with the medium and the interaction potential between the c-(anti)quarks (spurious). The modifications from elastic collisions do not depend on the size of the time step because elastic collisions are point-like in time. The variations of the total momentum of the c-(anti)quarks due to the interaction potential for a given time step will depend on the value of the time step.

The EPOSHQ code that executes our model is based on the trajectories of all c and c-quarks in Bjorken time, which serve as clock for the hydro evolution as well as for the HQ, so a shorter Bjorken time step makes the trajectories smoother. Once EPOSHQ generates the trajectories, the Minkowski time step used by our ensuing specific part of the code does not influence the results.

The cheapest way to explore the influence of the Bjorken time step value is through the behavior of the collisional rate when the elastic collisions with the medium are not active. This allows us to know that part of the total collisional rate which is spuriously created by the interaction potential only. Secondly, by modifying the Bjorken time step (as displayed on figure 5.12), we can study the effect of the time step and check if there is a convergence when ∆τ → 0. In figure 5.10, we have displayed the collision rate for different values of η cut and medium elastic collision and c c interaction potential. It should however be noticed that when we employ only the interaction potential (no elastics collision allowed), the absolute value of the collision rate is much lower (around 10%) compared to the full collision rate. This is due to the fact that in most of the time steps the momentum transfer due to the potential interaction is smaller that η cut (as far as η cut has a finite value) and therefore these time intervals do not contribute. This comes together with the fact that interaction potential can only make significant variations in the momentum of those c-(anti)quarks which are really close to pair partner.

We construct our trajectory of the heavy quarks in Minkowski space from their trajectories in the hydro-dynamical part of the EPOS calculation, which are given in Milne coordinates. Therefore a shorter Bjorken time step allows for a smoother interpolation which reduces subsequently the momentum transfer during a time step in Minkowski space. Therefore our results depend on the Bjorken time step, as can be seen from figure 5.12. Although the uncertainties derive from any particular selection of the Bjorken time step will be relatively small, of the order of

∼ Γ potential coll (∆τ) -Γ potential coll ( ∆τ 2 ) Γ full coll (∆τ)
Where Γ potential coll

is the time integrated collision rate with only interaction potential active (no medium elastic collision as the calculations in the figure 5.12), Γ f ull coll is the time integrated collision rate with both interaction potential and medium elastic collision (similar to the results of the figure 5.10) and ∆τ is a typical Bjorken time step (in our simulation ∆τ=0.1 fm/c).

On the other hand, the Minkowski time step, although subordinate in hierarchy to the Bjorken time steps, plays a role (in a similar way that Bjorken time) in the collision rate when only interaction potential is considered medium elastic collision are neglected.

Diminishing the Minkowski time step also diminishes the collision rate absolute value when there is only potential interaction, but with physical limitation of ∆t > ∆τ. Which comes from the relation ∆t = ∆τ 2 + z 2 which imposes the hierarchy of the Bjorken time step over the Minkowski one.

Role of the Q Q potential and/or the interactions with the medium Finally we compare the integral rate results in Fig. 5.13 switching ON and OFF the interaction potential and medium collision (ON and OFF). The trivial case, interaction potential and medium elastic collision being set OFF, yields a zero collision rate.

The maximum value ( of the integrated collision rate) is shifted towards large p T due to the absence of medium elastic collision. The comparison between the cases in which we always have medium elastic collision but potential is ON and OFF is more challenging. Not only because the fact that when interaction potential and medium collision are (ON) is lower (in absolute value) but also the significant different with respect to the case of potential (OFF).

To address this behavior, we consider two effects, first when potential is active, those pairs trapped in the potential well will be more difficult to dissociate (need more energetic gluon to break the binding process). Meaning that the evolution of those pairs will evolve toward larger possibility of charmonium formation as time goes on, which makes collision rate less "negative". On the other hand, when medium collisions take place, those pairs with very low p T (and mid-rapidity) come closer and closer to thermal equilibrium. These multiple collisions prevents those pairs to leave the high density zone . This can significantly helps that the c and c get trapped in the potential well.

Final production probability

After scrutinizing each of the three processes that contribute to the total probability of quarkonium formation we are ready to present their influence on the total probability under the previously discussed scenarios.

Our master equation for the total probability was already established in chapter 3 (A.4), and from the discussion in the previous sections of this chapter 3, we learned (in general words) that primordial production is highly concentrated in the earliest phase of plasma evolution (although not at a single point in time). This contribution is not much influenced the potential interaction between c and c due to the high screening at that stage. The time-integrated local rate feels strongly the interaction potential because it has contributions from later time in the evolution of the plasma when the temperature is lower. Both, primordial and local rate, are strongly modified by elastic collisions with medium particles, which shift their p T -spectra towards lower values. Finally, the collision rate, which is obviously modified by collisions with the medium, is, as the the local rate, sensitive to interaction potential. Unlike the two previous contributions, it also has a weak dependence on the Bjorken time step through the influence of the potential interaction on the total momentum of the c-(anti)quarks. In this section, we will study how all these features pile up together and the consequence for the production of the J/ψ states.

Time evolution of final production probability

The time evolution of the total direct production probability is a fascinating study that describes how the probability of total production evolves as a function of the evolution time of the QGP in the computational frame. In figure 5.14, the results of the temporal evolution of the probability of formation of J/ψ for different centrality bins is presented. The production probability shows an increase at early times due to the arrival 5.6. Final production probability 137 of more and more c and c quarks, which enter regions with T < T diss ) and hence where the formation of bound states (primordial production) is possible.This growth is because each "new" c or c available to form bound states adds to the number of c or c in free plasma and regions with a local temperature lower than T diss . This increasing behavior continues until all c and c, created in the initial hard collisions, arrive in regions with T < T diss ).

From this point on, history differs depending on whether the interaction potential is active or not. If there is no interaction potential, the dynamics favors dissociation over the recombination processes. This means that we lose a fraction of the primordial production. While for the case of an active interaction potential, the situation is the opposite. Here, after the primordial production, we observe that recombination dominates over dissociation. This is deeply related to those pairs trapped in their respective potential well, which often seem to resist getting dissociated by thermal kicks, especially at a later stage of the evolution when the interaction potential is stronger and the energy-range of those thermal kicks is narrower. The later should be interpreted as some interactions with a dissociative nature do not manage to dissociate the pair (separate the c from c in phase space) because the interaction potential is restored fast enough (before the c and c quarks move away from each other) and can compensate the changes produced by the dissociative processes.

It is also important to notice that none of the evolutions presented in figure 5.14 supports the adiabatic picture (which would correspond to a frozen evolution) which was implicitly advocated in early studies of J/ψ suppression. With the c c interaction ON, one observes in fact a steadily raising number of bound J/ψ, which signals that our coalescence mechanisms is reinforced with evolution time (as described above), contrarily to the instantaneous coalescence mechanism that is usually applied to describe parton hadronization.

p T -Spectrum of J/ψ Final Production Probability

Continuing our analysis, the result for the final direct production probability p T -spectra are presented for the different scenarios of interest. First, the centrality dependence of p T -spectra is shown in a similar way as in the previous sections. The pertinent feature of this set of results, displayed in figure 5.15, is the overall effect of the interaction potential on the final J/ψ production probability. The interaction potential increases substantially the formation of bound states by keeping the pairs with very small phase space relative coordinates (p rel , r rel ) together and by protecting them from dissociation by thermal gluons and light quarks. The increase of the production probability with the centrality is a direct consequence of the increase of the initial c c pairs produced in the initial hard scatterings due to the increasing number of participants. In figure 5.16, the case of the most central collisions is further scrutinized by comparing the respective roles of Q Q potential and interactions with the medium. ALICE data are also reported on this figure and are found to be in good agreement with our "full" model (both types of interactions turned ON).

In the cases in which medium elastic collision are ON, the production probability at very low p T (p T ≥ 2 GeV/c) is larger than the cases in which medium elastic collision are neglected (OFF). This is independent on whether the interaction potential is active or not. Another consequences of the medium elastic collision, the overall shifting of the distributions towards low p T , is still observed at the final stage of the evolution. When we turn the medium elastic collision OFF, c-(anti)quarks travel through the medium as if they were white objects (with the exception of the case of med coll (OFF) and potential interaction (ON) in which they can interact with their corresponding antiquark only), without undergoing collisions. Which implies that they move in a ballistic regime, similar to the case of pp collisions. And also analyzing the results of figures (5.8, 5.8, 5.13 and 5.16), we see that in the cases in which the medium elastic collision was set to OFF the final production probability is mostly defined by the primordial production due to this ballistic regime.

From this analysis, one may conclude that at low p T , the medium elastic collisions help the quarkonium state formation by "delaying" (keeping the c-(anti)quarks closer for longer time scales) the pairs separation process. While at intermediate and high p T medium elastic collision contribute to the suppression of the production by shifting the c-(anti)quarks population towards lower p T values due to the energy loss processes.

Another significant feature is that the interaction potential increases the total production of J/ψ by an approximate factors of 3 and 4 independent of whether the elastic medium collisions were ON and OFF. The interaction potential affects only c c pairs with small relative coordinates (especially in position space). These pairs are consequently the pairs that contribute the most to the total probability. The rest of the pairs (neither close in position space, or in momentum space, or in both) are exponentially suppressed.

The action of the interaction potential causes the c and c of a given pair to get closer and closer as time advances or at least extends the time that c and c remain close. This increases the probability exponentially . Although the number of pairs in which the potential action is significant, is really small (only the pairs in which the relative distance is of the order of the inverse of the Debye mass r rel ≈ 1 m D ), they are the most important contributors to the overall production probability.

J/ψ nuclear modification factor

One of the main objectives of the thesis is to present our predictions for the nuclear modification factor R AA .

R PbPb = d 2 N Final PbPb dp 2 T d 2 N pp-scaled dp 2 T .
(5.9)

In our approach the "pp-scaled Monte Carlo" production will be used for the denominator. Its calculation has been explained in section 1 of this chapter.

To begin with, the individual contributions of each of the three terms that make up the final probability of production of J/ψ for the set of 4 scenarios are shown in figures 5.17 and 5.18. They are obtained from medium elastic collision and potential interaction switched (ON) and (OFF). To address each individual contribution properly, the concept of partial nuclear modification coefficient needs to be introduced. For example, the primordial partial R AA would be the application of the equation (5.9) considering only the primordial production As all contributions are rescaled by the same scaled pp production, they can be directly inter-compared. If one switches off the elastic collisions with the medium (OFF), the contributions of both rates (collision and local rate) are extremely small or directly null (in the case of the collision rate when also interaction potential is off (OFF). As anticipated in the previous section, the primordial production is then very close to the final production. When the elastic medium collisions are switched on (ON), the situation is more complicated due to the more active role of both rates. As already mentioned, the integrated contribution of both rates is of the same order and when there is no potential interaction the suppression is stronger. To conclude the section, the comparison with the experimental data reported by the ALICE collaboration for the mid-rapidity nuclear modification factor for PbPb collisions are displayed in Fig. 5.19 for two different ranges of centrality (0 -20% and 20 -40%) For low p T our predictions are in the vicinity of the experimental data and previous model predictions [START_REF]Centrality and transverse momentum dependence of inclusive J/ψ production at mid-rapidity in Pb-Pb collisions at √ s = 5.02 TeV[END_REF]. For intermediate and large p T our calculation yield larger values than the experimental data.

Two major factors may be responsible for this disagreement. First, our model for pp-scaled production clearly underestimates the experimental data at large p T values, with direct consequence for the nuclear modification factor. Besides, the extra contribution one should include at large p T (c c from gluon splitting) would undergo extra energy loss as compared to the one considered up to now in the model, which could lead to extra suppression.

Second, when the interaction potential is active and we include both, collision and local integrated rates, the behavior of the p T -spectra at high p T slightly overshoots the experimental data. This effect is directly related to the action of the potential for that sector of p T . The description of the relativistic bound state formation has proven to be a very difficult task, mainly due to the challenge to find a way to map the many-body relativistic dynamics from the laboratory frame to the specific c c center of mass. The adopted algorithm supplemented by the Rothkopf potential, is expected to be very good at low p T , but its performance at intermediate and high p T should be further improved, while the question of the c c potential at high p T is a hot topic under discussion in the HIC community. Finally, the behavior of the R AA as a function of the average number of participants in the PbPb collisions is presented in figure 5.20. This result allows us to extend the study of the R AA by the different bins of centrality. Therefore it allows us to see the behavior of our model while important initial parameters vary, such as the initial number of c c pairs, the geometry of the collision, and the initial conditions of the plasma. As previously commented in this chapter, the results presented for J/ψ production correspond, in our model, to direct production and are not affected by any CNM effects. Therefore, when comparing results with the experimental data from ALICE collaboration corresponding to inclusive production, some underestimation should be expected, as we can see for the cases of the semi-central collision in figure 5.20. Even so, in cases of extreme centrality (both high and low centrality), our predictions overestimate the experimental results, especially at high p T (see figure 5.16) where a further contributions should be considered due to the decay of excited charmonium and B-meson (which become more relevant at high p T ).

The c c pairs in our calculations are produced in an uncorrelated way, so the high p T production is conditioned by the statistical probability of finding the highest collinear c c pairs that the interaction potential can trap. This behavior is purely statistical, which may explain the fluctuating nature of our predictions (overestimation for extreme centrality cases and underestimation for semi-central cases). The high p T sector is especially sensitive to this statistical fluctuation since it naturally has less 5.8. Diagonal and off-diagonal contributions 145 population than the low p T sector. So the impact of fluctuations will be relatively more significant in the high p T sector than the low p T . From a technical point of view, running an even larger number of simulations (increasing the statistics) can cure these fluctuations. But it is the physical setting where the most exciting aspects are found. Comparing the results presented in this thesis (uncorrelated c c production) with corresponding results with correlated c c pairs and considering cc pairs coming from jets and gluon splitting can generate an interesting debate on how the production of charmonium is fed in the different regions of p T .

All together, the comparison presented in figure 5.20 should be taken with a grain of salt due to the uncertainties that are still affecting the c c production cross section in pp as well as the CNM effects.

Diagonal and off-diagonal contributions

One of the most important advantages of our model is that for each pair at each instant in time the probability of forming a J/ψ (seen as a bound state of the local Hamiltonian) can be calculated. This allows, among other things, to study how those pairs contribute individually, in particularly those formed by a c and c from the same production vertex. We may think of a matrix arrangement where indices i and j run over the production vertices at which the c(i) and the c(j) are produced. The pairs from the same vertex will have the form (ii) and we called their contribution "diagonal contributions" (such "diagonal" contribution corresponds to what other models name "primordial component"), while the rest of the form (ij) can be referred to as "off-diagonal". The relationship between both contributions (which are consistently described within our approach) is of great interest because it directly describes how the final probability is composed. Analyzing the results shown in figures A.3 and 5.22, we can see that the diagonal contribution increases with increasing p T . This can be explained by the combination of two effects. On the one hand, the contributions to the final total probability at high p T are concentrated in those extremely colinear pairs (pairs with very small relative momentum) and relatively close in position space. Although the initial distribution of the produced c-(anti)quarks is uncorrelated, there is some statistical probability of forming charmonium states from a highly colinear pair and with a small relative distance. These highly colinear and close pairs in position space will strongly feel the potential for c c interaction. The latter is because although the energy of the center of mass in the computational frame is very high (relativistic), the interaction potential works in the center of mass frame, where for the case of interest (very colinear and close pairs), the dynamics are non-relativistic. In this scenario, the action of the interaction potential is substantial, keeping c and c close for a long time interval. In the low p T sector most c and c quarks have a high degree of thermalization, so they do not preserve the memory of the origin of their production (they are independent of the vertex of origin). This means that the contribution fraction from the same vertex (diagonal) is purely statistical, especially for later times when the degree of thermalization becomes larger.

A practical way to study how both contributions (diagonal and off-diagonal) participates in the total probably is by introducing some Ratio as indicated in the equation (5.11)

Ratio = 1 + d 2 N off diag dp 2 T d 2 N diag dp 2 T (5.11)
The results obtained in the study of this auxiliary quantity are presented in the figure 5.23, where two different scenarios have also been prepared, marked by the labels: X1 and X2. The X1 label refers to the fact that during the simulations carried out by EPOSHQ, the c c cross-section production considered corresponds to that reported by the different experimental collaborations, which yields to dN c c dy = 15.6656 at 0-20% centrality interval. Whereas the X2 label refers to the cross-section being duplicated, which then will yield to dN c c dy = 31.6373 as indicated in the table 5.1, which correspond to our standard choice for all the results presented in this chapter.

The cause of the decrease in both ratios when passing from X2 → X1 was already mentioned above, but it must also be taken into account that our interaction potential more strongly modifies the high p T production. , taking inspiration from the concept of nuclear modification factor, which was already discussed in the previous section.

R diag/off-diag AA = d 2 N diag/off-diag PbPb dp 2 T d 2 N pp-scaled dp 2
T .

(5.12)

As we can see in the figure 5.24, the off-diagonal contributions completely dominate the production at low p T and the c-(anti)quarks have a high degree of thermalization. So pairs from off-diagonal combinations (c and c from different vertex) are statistically more likely.

In the same way, the production at high p T is exclusively concentrated in highly collinear pairs coming from the diagonal contributions (implying that they are relatively close), which are the pairs most susceptible to the interaction potential.

To conclude this section, the results corresponding to the nuclear modification factor of the diagonal and off-diagonal components (5.12) of the final production at mid-centralities (20-40% and 30-50%) are presented in figure 5.25 . These results will serve as a reference for the discussion of the results of the diagonal component of elliptic flow (v2) in the following section and complement the results presented in figure 5.24. Analyzing the results of figure 5.25, we can verify (as mentioned before) that as the production of c c pairs decreases, the diagonal component of the production becomes more relevant. This behavior will continue to the limit in small systems such as pp collision, in which the total probability is almost entirely made of diagonal contributions.

J/ψ Elliptic Flow

The study of the flow coefficients defined at eq. (2.38), especially of the elliptical flow v 2 , is of great interest, because they are expected to be contain information about the quarkonia production times. If they are produced early, before the heavy quarks acquire v 2 in the collisions with light quarks and gluons, the v 2 of the quarkonia would remain small because the interaction of a color neutral object with the medium is weak.

We study here how each of the three terms that make up the final probability depends on one or more medium properties, such as the local temperature and the influence of collisions with bulk particles. Our model is capable of calculating at each given time the impact that any modification or interaction, suffered by c and c quarks, has on the total probability and its repercussion with respect to the rest of potential partners for the formation of pairs. This guarantees that any collective effect that affects several c and c quarks is reflected in the final probability of production. Therefore, we consider that the greatest strength of our model may lie in the predictions of these collective phenomena.

In figures 5.26 and A.4, the experimental results for the elliptic flow, reported by the ALICE collaboration [START_REF]J/ψ elliptic and triangular flow in Pb-Pb collisions at √ s =5.02 TeV,ALICE Collaboration[END_REF] are compared with our predictions for the case in which the potential interaction and the medium elastic collision were set ON. From this comparisons we can conclude that our predictions agree relatively well with the ALICE results . We see in figures A.4 and 5.26, that our model of a dynamical production of producing composite objects (charmonium) gives reasonable results while considering only the direct component of the J/ψ production and under uncorrelated regime of c c. In this scenario, the model offers better results than most of the models available on the market [START_REF]φ and Ω production in relativistic heavy-ion collisions in a dynamical quark coalescence model[END_REF][START_REF] Molnar | Elliptic Flow at Large Transverse Momenta from Quark Coalescence[END_REF] by providing a good agreement both with the nuclear modification factor and the elliptic flow. The direct dependence on the number of components available for forming bound states and their respective quantum states is the perfect ingredient to study collective effects under the HIC environment. Unlike other models that deal with recombination processes such as SHM, which assumes thermal and chemical equilibrium, our model only relies on the distribution in the phase space of the c and c quarks, which decides about the probability of quarkonium formation based on the projection of the two body Wigner density on the quantum mechanical ground state of the quarkonium.

Taking advantage of the virtues of our model and starting from the definitions made in the previous sections where the differences between the diagonal and off-diagonal contributions were discussed, we present, in the following figures, the corresponding results to the elliptic flow of the diagonal component of total production. From the results presented in the figure 5.28 we can conclude that the diagonal component of the total production contributes substantially to the total v 2 while being smaller then the off-diagonal v 2 . We recall that in our model, c and c are produced in an uncorrelated way2 and they interact as independent colored objects when they propagate in the hot medium. This scenario also includes the diagonal pairs, which are allowed to interact with the medium particles at any stage of the evolution (until they abandon the plasma forming any bound state D-mesons, J/ψ, ...) inside the QGP even if there are relatively close to a potential partner.

This approach arguably offers a different vision than other models on the market, which consider that these diagonal pairs are practically white objects, what completely cancels the interaction of these pairs with the medium's particles. Without these interactions, the said pairs will not be able to absorb the flow of the medium, what will therefore result in zero or very little3 contribution to the total v 2 from the diagonal component of the production.

From a more theoretical viewpoint, it is known (see for instance [START_REF] Blaizot | Quantum and classical dynamics of heavy quarks in a quark-gluon plasma[END_REF]) that compact white objects do not interact with the QGP, while the interaction rate of those objects with the QGP increase quadratically with their size r sing , until r sing ∼ l corr , the correlation length of gluon thermal fields. From this value on, both Q and Q interact independently. In this respect, one may consider that the traditional models based on a non-interacting primordial singlet component and our approach explore both facets of the more involved reality, which leads to a bracketing of the experimental elliptic flow.

Chapter 6

Conclusions and perspectives

Conclusions

In this thesis we have studied the question whether the attractive interaction between a c and a c can yield an enhancement of the R AA at LHC as compared to RHIC energies. This has been experimentally discovered and the question has been raised whether the much higher number of c c pair, created at LHC, may be at the origin of this enhancement

The starting point of our work is the model of Remler who proposed to calculate the pair production by studying the probability that a c and a c quark form a quarkonium. This probability is calculated before and after each collision of a heavy quark with the constituents of the plasma. The difference between both probabilities increases or lowers the probability that a quarkonium is created in the QGP.

Adapting this low energy model to the high energy scenario of ultrarelativistic heavy ion collisions several problems had to be overcome: the numerical implementation of a multi-particle time evolution scheme, the introduction of a time dependent density operator (Wigner function) referred to in this works as a "local basis" correction, Q Q potential interaction, and corrections due to relativistic dynamics.

The tasks completed in this project show that it is possible to develop a semi-classical model derived directly from the most fundamental principles of quantum mechanics that allows to observe the temporal evolution of the probability of the formation of bound states without assuming particular conditions on the surrounding environment. This was possible thanks to the solid basis on which the [109] model was built and the great success it has had in the different scenarios in which it has been tested [109,122].

This work presents a first approach to the production of quarkonia. The model can Chapter 6. Conclusions and perspectives be extended in several ways. The relativistic kinematics of the potential interactions is employed in first order in β where time retardation effects are not present. The trajectories of the c and c are calculated classically and the transition between Milne and standard coordinates can be improved. The results show, however, that a potential interaction between c c pair enhances R AA in the region where it is observed experimentally and therefore our approach demonstrates that this interaction may be well at the origin of the different observations at RHIC and LHC. Our model also qualitatively describes the experimentally observed v 2 of J/ψ.

Perspectives

As mentioned in chapters 3 and 5 of the thesis, the interaction potential used in the simulations is the Rothkopf-like potential, which has shown excellent success for low energies (low p T ). Still, its application for high p T remains to be consolidated. Another problem related to relativity is the retardation effects. These effects appear when the relative momentum is high, thus causing the action of the interaction potential to be retarded. In our calculations, the interaction potential acts instantaneously (in the center of mass frame), meaning that we intrinsically assume that although the total center of mass energy of the pair (quarkonium) is relativistic, the relative momentum of the pair will be small (not relativistic). This approach is very realistic in our model because, given the mathematical nature of the Wigner density function, those pairs with high p rel will be exponentially suppressed in the probability counting. Yet, the extension of the retardation effects was not explored here and may also be affecting our prediction.

During the realization of this project, the possibility of implementing a electromagnetic potential (chromo-electromagnetic like version ) with a retardation correction of the Lienard-Wiechert type was explored, but no fruitful results were achieved 1 , so this idea was discarded at that time by a more solid and practical option (Rothkopf potential). But this idea may be pursued. Our model can be a basis for examining such ideas and their influence on observables such as R AA or v 2 .

Another aspect that can be revisited is the influence of the local medium temperature in the Wigner density function (time-dependence Wigner function). In the present 1. Due to the non conservation of the Q Q energy measured in the center of mass and the associated collapse of the state. work, the time dependence results in a temperature dependence with enters in the Gaussian width σ = σ(T) = ( 23 < r 2 (T) >) 1/2 . The numerical dependence was derived from the behavior of the mean square radius of the J/ψ at ground state under a Rothkopf interaction potential.

A more rigorous way to derive this temperature dependence may be obtained by deriving the exact solution of the Schrodinger problem under a given temperature potential (like Rothkopf potential) and by obtaining the corresponding Wigner density. In case one wants to preserve the double Gaussian approach for a Wigner density function, a more sophisticated temperature dependence may be considered. It is also important to emphasize that from the perspective of the interaction of the HQs with the QGP (for example, medium elastic collision), variables such as relative distance should be taken into account when processing the interaction. This analysis will allow distinguishing between pairs whose relative distance is very small (and therefore, their interaction rate should be relatively low) from those whose relative distance is large (in which case the interaction rate should be relatively higher).

Finally, another perspective to add is that even in a "crude" state, our model can be used to study the formation of different excited quarkonium states as well as the formation other heavy flavor state formation like recently observed B + c mesons [133]. Of course, several modifications need to be implemented to carry out the desired simulations: Fits to a Wigner density function (including the necessary quantum number correction) corresponding to the desired excited state, its dissociation temperature, and the mass of the excited state. This is, in turn, one of the most interesting perspective to pursue in the near future. conditions are fulfilled. We can make now a second variable transformation P 0 = √ s CoshY, P = √ s SinhY, q 0 = QSinhy, q = QCoshy (A.14) with d 2 Pd 2 q = √ sd √ sdYQdQdy, including both, the q = +QCoshy and the q = -QCoshy sector. We assume now that f i depends only of the invariant relative momentum Q = q 2 -(q 0 ) 2 f i (Λ -1 (q 0 , q)) = f i (Q).

(A.15)

With

δ(P • q) = δ( √ sQSinh(y B ± Y)) = δ(y B ± Y) √ sQCosh(y B ± Y)
where the ± stand for the q = ±QCoshy sectors, respectively. We can now define our relativistic state

|Ψ y B ,i = dQ √ s dYδ( √ sSinh(y B -Y))( f i (+Q) + f i (-Q))
( f i (+Q) + f i (-Q))( f i (+Q) + f i (-Q)) = dκdY √ s δ(y B -Y)δ(y B -Y) sCosh(y B -Y)Cosh(y B -Y) f i (κ) f i (κ)× = δ(y B -y B ) dκ s 3/2 (κ) f i (κ) f i (κ) (A.18)
and we can derive the invariant orthogonality relations provided we require

dκ s 3/2 (κ) f i (κ) f j (κ) = δ ij .
By introducing the non relativistic wave function

Φ i (q) = f i (q) s 3 4 (q)
(A. [START_REF] Balk | Heavy Quark Effective Theory at Large Orders in 1/m[END_REF] we obtain the non relativistic orthogonality relation dqΦ(q)Φ(q) = δ ij .

Although we have been able to develop a recipe for the construction if an orthogonal invariant base, there are still some comments to be made: first, we do not strictly recover the desire normalization Φ y B ,i |Φ y B ,i = 2Eδ(pp )δ ij . This can be explained due to the fact that for a given state i and a given rapidity y B , the total energy depends on the relative momentum κ. This has consequences for the orthogonality relation A. [START_REF] Balk | Heavy Quark Effective Theory at Large Orders in 1/m[END_REF], which slightly differs from A.6. However, the relationship A. [START_REF] Golden | On the Hilbert space of the heavy quark effective theory[END_REF]

Résumé de la thèse en Français

Introduction

Ces dernières années, la recherche dans le domaine des collisions d'ions lourds (HIC) a été caractérisée par un travail intense dans le domaine expérimental, réalisé par les collaborations associées au collisionneur LHC (ALICE, CMS, LHCb et ATLAS), ainsi que les collaborations associées au collisionneur RHIC (PHENIX et STAR). Dans le domaine théorique et phénoménologique, les discussions sur les différents modèles et les approches qui tentent d'expliquer les processus propres aux HIC continuent d'être très débattues, malgré les nombreuses données expérimentales que les différentes collaborations ont collectées. En particulier, les débats sur les modèles se référant à la production des différentes particules n'est de loin pas terminée. Malgré les discussions animées, un certain consensus a été atteint concernant quelles observables sont les plus adaptées pour explorer la dynamique des HIC. Parmi meilleures observables, la production de particules s'impose comme l'une des observables les plus prometteuses, en particulier celles qui proviennent des phases initiales des HIC. Par conséquent, créer des modèles théoriques et phénoménologiques capables de prédire la production de particules dans les HIC est une priorité de recherche pour la communauté des collisions d'ions lourds. L'objectif principal de ce travail est de présenter et de discuter un modèle de coalescence dans l'espace des phases qui permet de prédire et d'explorer la production de quarkonia lourds (états liés d'une production de quark lourd Q et de son antiquark correspondant Q ¯ ) dans les HIC ainsi que dans les collisions proton-proton. Notre intérêt particulier pour les quarkonia lourds est basé sur le fait qu'une fraction importante des quarkonia est produite dans la phase précoce d'une collision d'ions lourds, puis traverse le plasma quark-gluon (QGP), et quitte finalement le plasma sans atteindre la thermalisation complète (donc ils ne perdent pas tout souvenir des processus auxquels ils ont participé). En étudiant la formation de quarkonium avec notre modèle, nous pouvons obtenir des informations sur la façon dont les différentes interactions des systèmes de quarkonium (Q Q) avec les particules du milieu modifient la probabilité de production, y compris l'effet des variables QGP telles que la température et les caractéristiques macroscopiques du plasma telles que les effets collectifs (de type flot).

Formalisme de Remler

Les idées de notre modèle ont été postulées par E. A. Remler [109]. Il a résolu les équations d'évolution pour l'opérateur de densité quantique en utilisant les densités de Wigner. Ceci permet de calculer le taux de production de particules composites. Ce formalisme a ensuite été appliqué avec succès au problème de la formation de deutérium dans les réactions nucléaires [110]. Plus récemment, une approche de coalescence s'inspirant du modèle de Remler (mais n'adressant pas l'aspect de l'évolution temporelle) a été appliquée avec succès aux production de charmonium dans les collisions proton-proton [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF]. Notre objectif principal est de construire un modèle basé sur le formalisme de Remler et applicable aux conditions créées dans les HIC. Toutes les informations sur un système peuvent être extraites de l'opérateur de densité à N corps du système (ρ N (t)). L'opérateur densité du système complet obéit à l'équation de von Neumann [110] :

∂ρ N /∂t = -i[H, ρ N ]
Où H est l'hamiltonien du système complet et peut être décomposé comme suit :

H = Σ i K i + Σ i>j V ij
où K i est l'opérateur d'énergie cinétique de la particule i et V ij est l'interaction entre les particules i et j. Les mésons, comme un J/ψ, sont des objets à deux corps décrits par l'opérateur de densité ρ Φ = |Φ >< Φ|. où Φ est l'état quantique correspondant à les mésons lourds formé d'un Q et d'un Q. Ainsi

P Φ (t) = Tr[ρ Φ ρ N (t)],
où la trace est prise sur toutes les coordonnées N-corps (qui incluent les degrés de liberté de Q et de Q), mesure la probabilité de trouver le Q et le Q dans un sous-espace défini par ρ Φ à l'instant t. Il est difficile d'utiliser l'éq. (3) directement dans certaines implémentations numériques de type "cascade". En effet, dans ce type d'approche, la distance entre les quarks Q et Q à la fin de l'expansion QGP est grande par rapport au rayon du méson Φ et donc P Φ tend vers zéro. Pour surpasser ce problème, nous calculons le taux de formation et de décroissance de Φ, Γ Φ (t), qui est liée à la probabilité de trouver un J/ψ au temps T, P(T), en réécrivant l'équation (3) comme La fonction d'onde relative initiale ou densité de Wigner de la paire Q Q, peut être calculée à partir des premiers principes en résolvant le problème de Schrödinger exact pour un potentiel. Ainsi en s'inspirant de l'approche de l'oscillateur harmonique, on peut ensuite approximer la fonction de Wigner du système à deux corps par une double gaussienne [START_REF]The production of primordial J/Ψ in p+p and relativistic heavy-ion collisions,Taesoo Song,Joerg Aichelin and Elena Bratkovskaya[END_REF] par rapport à la distance relative r et à la quantité de mouvement relative p de la paire dans le centre de masse: Effets de la température sur la production de quarkonium, la base locale En tenant compte maintenant que les quarks lourds sont dans un milieu de température T, l'hamiltonien H 1,2 dépendra désormais également de la température. Ceci induit que ρ Φ et donc W Φ dépendront également de la température, puisque la relation [ρ Φ (T), H 1,2 (T)] = 0 doit être conservée.

P Φ (T) =

Le taux local

En revenant dans notre analyse à l'équation du taux et en tenant compte maintenant que ρ Φ (T) dépend du temps par la température T = T(t). En conséquence de cette considération, l'équation d'évolution et le taux présente maintenant un deuxième terme Ainsi, on peut soutenir que l'élargissement de la dimension de l'état lié (rayon carré moyen) avec T est responsable de ce comportement. Et c'est le chemin par lequel la température modifie la base de l'état lié (base locale) dans notre approche.

Γ eff = Tr[ρ Φ (r,

Production initiale de quarkonium

Un autre effet de la température est que si la température locale du plasma est supérieure à la température de dissociation de l'état quarkonium, ledit état quarkonium ne peut pas être produit ou s'il existe déjà il sera détruit par l'interaction avec les particules du quarkonium avec le milieu. Par conséquent, dans notre approximation, la production initiale de quarkonium sera distribuée dans le temps, car nous devons attendre que les quarks lourds arrivent dans les régions du plasma moins chaude, où la formation d'états liés est possible. Comme la distribution de température dans le plasma est anisotrope, les quarks lourds seront prêts pour la production de quarkonium à des instants différents t Q, Q init .

Cela dit, la probabilité d'une production d'état de quarkonium à partir d'une paire Q Q aura la forme

P Q Q(t) = P prim (t Q, Q init ) + t t Q, Q init (Γ coll,Q Q(t ) + Γ loc,Q Q(t ))dt
Notons que cette équation ne fait référence qu'à la probabilité de formation d'un état quarkonium à partir d'une seule paire Q Q. Pour connaître la probabilité totale, il faut additionner l'équation précédente sur toutes les combinaisons possibles de Q et Q présentes dans le plasma à l'instant t.

Extension relativiste des fonctions de Wigner et interaction des paires Q Q dans le plasma

Bien que le programme de Remler ait été mis en oeuvre à l'origine en association avec certaines évolutions de type "cascade" de la densité de N corps prise au niveau semi-classique, la véritable l'évolution semi-classique de W N s'avère avoir un rôle certain lorsque les particules sont sous l'influence du milieu pour des temps plus grands que le temps de Heisenberg (qui est le temps typique pendant lequel Q et Q se cyclent l'un l'autre dans une résonance donnée). Ainsi, l'algorithme de Remler n'est qu'une moitié du modèle, mais il doit être complété par l'autre moitié qui consiste en une évolution semi-classique réaliste des quarks lourds au sein du QGP comprenant une interaction à la pointe de l'état de l'art entre Q et Q avec une extension relativiste de la fonction de Wigner.

Extension relativiste de la fonction de Wigner

La formulation relativiste des états liés du quarkonium est une tâche difficile, car nous ne pouvons plus séparer le centre de masse et le mouvement relatif. En outre, un effet de retard peut apparaître si la quantité de mouvement relative de la paire est grande dans le centre de masse, tandis qu'en plus les temps dépendent du système de référence. Au cours du développement de notre travail, une extension relativiste de la fonction de Wigner a été dérivée, qui s'écrit sous la forme W i,u lab T,Φ ,y Φ (y, u T ; r lab , p lab ) = 1 (2π) 3 δ(yy Φ )δ 2 (u T,Φu T )× W i,NR (r cm (r lab ), p cm (p lab ))

où W i,NR , indique la version non relativiste de la fonction de Wigner, qui sera évaluée dans le centre de masse, alors que u T,Φ et y Φ sont la composante transversale de la vitesse et la rapidité du centre de masse par rapport au référentiel du laboratoire (lab).

Évolution relativiste des paires Q Q

Chaque état Q Q sera donc propagé en résolvant le problème à deux corps dans le centre de masse (cm) puis en effectuant des allers-retours systématiques de type "transformation de Lorentz" entre le référentiel du centre de masse (cm) et celui du cadre de calcul global (laboratoire). Ceci constitue une autre différence entre l'algorithme actuel et le formalisme original de Remler, qui est évalué dans un cadre unique (laboratoire). La solution du problème à deux corps a été réalisée dans le cadre de notre correction de l'algorithme d'origine. C'est parce que lorsque nous avons affaire à un système relativiste Q Q, deux problèmes majeurs dérivés de la conservation de l'énergie et de la quantité de mouvement voient le jour. La non-conservation de l'énergie lors de l'évolution d'un Q Q dans un référentiel spécifique (par exemple le centre de masse) peut conduire à des comportements tachyoniques, qui bien sûr peuvent détruire tout effort de calcul. D'autre part, des paires avec des énergies ultrarelativistes peuvent aussi comporter de grandes quantités de mouvements relatives (même dans le centre de masse), ce qui implique que certains des effets d'interaction retardées (interactions de type Lienard-Wiechert) sont attendus.

In fine, les équations qui décrivent le mouvement des paires sont exprimées sous la forme: mγr θ = L 2 mγr 3 qui décrit l'évolution de la variable angulaire θ, tandis que l'évolution de la variable radiale r est décrite par

L 2 d 2 dθ 2 1 r -γ C m + L 2 r = 0
où est γ le facteur de Lorentz du quark lourd mesuré dans le référentiel du centre de masse.

Intégration du modèle à la génération d'événements EPOSHQ

Malgré le fait que dans la section précédente les équations du mouvement (en particulier l'équation de la variable radiale) aient été obtenues en supposant un potentiel de Coulomb, dans nos simulations réelles, l'interaction forte entre Q et Q est mise en oeuvre via un potentiel évalué par le biais de calculs sur réseau par A Rothkopf, qui s'avère plus approprié pour décrire la dynamique des paires dans un milieu à température finie. La production, l'initialisation et l'évolution hydrodynamique des quarks lourds sont quant à elles réalisées par le générateur d'événements EPOS-HQ, générant ainsi les trajectoires classiques de ces quarks.

Notre modèle est exécuté en parallèle au sein de l'infrastructure EPOS-HQ. Une fois les quarks lourds produits et initialisés, nous suivons l'évolution de chacun d'entre eux au cours de leur trajectoire à l'intérieur du plasma. A chaque pas de temps, la méthodologie décrite dans les sections 3 et 4 est appliquée. Ce faisant, notre modèle est capable de prédire la production de quarkonium en fonction du temps, de la rapidity et de l'impulsion transverse pT, comme on peut le voir dans les résultats présentés sur la

Conclusions

Le point de départ de notre travail est le modèle de Remler qui a été proposé pour calculer la production de paires en étudiant la probabilité qu'un quark c et un c forment un quarkonium. Cette probabilité est calculée avant et après chaque collision d'un quark avec les constituants du plasma. La différence entre les deux probabilités augmente ou diminue la probabilité qu'un quarkonium soit créé dans le QGP. Adapter ce modèle basse énergie au scénario haute énergie des collisions d'ions lourds ultra-relativistes pose plusieurs problèmes qui ont dû être surmontés: la mise en oeuvre numérique d'un schéma d'évolution temporelle multi-particulaire, l'introduction d'un l'opérateur de densité (fonction Wigner) auquel il est fait référence dans ce document fonctionne comme une correction de "base locale", l'interaction potentielle Q Q, et des corrections dues à la dynamique relativiste. Les tâches réalisées dans ce projet montrent qu'il est possible de développer un modèle dérivé directement des principes les plus fondamentaux de la mécanique quantique qui permet d'observer l'évolution temporelle de la probabilité de formation des états liés sans supposer de conditions particulières sur le milieu environnant. Les résultats montrent cependant qu'une interaction potentielle entre les quarks c c améliore le R AA dans la région où il est observé expérimentalement et donc notre approche démontre que cette interaction pourrait bien être à l'origine des différentes observations au RHIC et au LHC. Notre modèle décrit aussi qualitativement les observations expérimentales du flot elliptique (v 2 ) des J/ψ.

FIGURE 1 . 1 -

 11 FIGURE 1.1 -Transverse momentum spectra of inclusive J/ψ (left) and ψ(2S) (right) mesons produced in pp collisions at √ s = 7 TeV. The blue triangles are data points coming from LHCb [2] collaboration and green triangles and red circles are data points from ALICE collaboration

FIGURE 1 . 2 -

 12 FIGURE 1.2 -Feynman diagrams of the HQ primordial (initial) production, light quark-antiquarks annihilation (bottom right) and the three leading order contributing processes of gluon fusion

FIGURE 1 . 3 -

 13 FIGURE 1.3 -Representation of J/Ψ production in pp collision according to the CSM .

FIGURE 1 . 6 -

 16 FIGURE 1.6 -ALICE inclusive J/ψ p T -differential cross section [20] at forward rapidity compared to mixed model NRQCD+FONLL(gray line) and NRQCD+CGC+FONLL (blue line).

FIGURE 1 . 8 -

 18 FIGURE 1.8 -Comparison between the QCD phase diagram (left) and water state diagram (right) [26].

FIGURE 1 . 9 -

 19 FIGURE 1.9 -Time scale of the different stages of AA collision

FIGURE 1 . 12 -

 112 FIGURE 1.12 -In the left panel, the nuclear modification factor as function of p T of averaged prompt D mesons in the 0 -20% most central Pb-Pb collision at √ s NN =2.76 TeV. Meanwhile in the right panel v 2 as a function p T of D mesons in 30 -50% centrality Pb-Pb collision at √ s NN =2.76 TeV. The results [44] are compared to model calculations implementing collisional energy loss.

FIGURE 2 . 2 -

 22 FIGURE 2.2 -Representation of the different process that take place in the heavy flavour in the QGP.

FIGURE 2 . 3 -

 23 FIGURE 2.3 -Partonic processes for the production of photons from (a) quark gluon Compton scattering, quark anti-quark annihilation process and (b) quark fragmentation.

FIGURE 2 . 4 -

 24 FIGURE 2.4 -The nuclear modification factor for inclusive J/Ψ production at √ s NN = 8.16 TeV [67] as a function of transverse momentum (Left) and as a function of the rapidity [68] for ψ(2S) compared to the corresponding quantity for J/ψ (Right).

FIGURE 2 . 5 -

 25 FIGURE 2.5 -Schematic representation of the Cronin effect.

FIGURE 2 .

 2 FIGURE 2.6 -R A F2 Representation of the different regime (effects) as function of the x values [78]

FIGURE 2 . 7 -

 27 FIGURE 2.7 -Mechanisms contributing to matter [82] induced changes in the yield of quarkonia: Debye mass effect (upper left), ionization by thermal gluons (lower left) and recombination (right).

FIGURE 2 . 8 -

 28 FIGURE 2.8 -Debye mass behavior with the thermal bath temperature [85].

FIGURE 2 . 9 -

 29 FIGURE 2.9 -Left: Debye radius representation in a cold medium ( relatively small Debye radius). Right: same for a hot medium (relatively large Debye radius).

FIGURE 2 .

 2 FIGURE 2.10 -J/ψ relatively contribution from its different state as function of the energy density

FIGURE 2 . 11 -

 211 FIGURE 2.11 -Inverse dissociation time (dissociation rate) [92] behavior with time.

FIGURE 2 .

 2 FIGURE 2.12 -Charm (left) and beauty (right) quarks friction coefficient comparison between the equation (2.29) (dashed lines) and Einsten's relation η D = κ L /2ET (continue lines) in the QGP[START_REF]Heavy-flavour and quarkonium production in the LHC era: from proton-proton to heavy-ion collisions[END_REF].In the results presented in the figure, the hard thermal loop (HTL) approximation was used to evaluated the soft interactions while the pQCD approach was employed to evaluates hard interactions with medium particles

FIGURE 2 . 13 -

 213 FIGURE 2.13 -Ratio of the number of J/ψ over the number of primordial produced ones (full curve) in a two-fireball [100, 101] model around mid-rapidity. The dashed (resp. dash-dotted) curve corresponds to the J/ψ yield from statistical coalescence at hadronization (resp. direct production with nuclear absorption and QGP suppression).

FIGURE 2 .

 2 FIGURE 2.14 -Prompt J/ψ yield in pp measured by LHCb at √ s = 7 GeV ([102]) compared to different theory predictions referred to as prompt "NLO NRQCD" (blue), "direct NLO CS" (gray), "Direct NNLO CS" (red) and Prompt "NLO CEM"(blue)

FIGURE 2 .

 2 FIGURE 2.15 -Typical sources of prompt J/ψ at low p T and mid p T ( [44]).

FIGURE 2 .

 2 FIGURE 2.16 -Typical sources of Υ(nS) at low p T and high p T

FIGURE 2 .

 2 FIGURE 2.17 -Differential cross section p T -spectra at √ s= 13 TeV for inclusive J/ψ experimental data from ALICE ([104]). In the left panel the models prompt and non-prompt J/ψ contributions are separated, while in the right panel the sum (inclusive) of both contribution are presented.

FIGURE 2 .

 2 FIGURE 2.18 -Representation of a non-central collision. The spatial anisotropy with respect to the x-z plane (reaction plane) leads to a momentum anisitropic distribution of the particle [105] in the final state

FIGURE 2 .

 2 FIGURE 2.19 -Time evolution of the energy density profile in coordinate space[START_REF] Snellings | Elliptic Flow: A Brief Review[END_REF] in a non-central collision. The z-axis is along the beam direction and the impact parameter is the same direction than x-axis.

FIGURE 2 .

 2 FIGURE 2.20 -p T -spectra D-mesons (black dots) elliptic flow compared with those of light hadrons reported by ALICE collaboration [79] in the left panel. While in the right panel the inclusive J/ψ v 2 p T -spectra at forward and mid-rapidity for PbPb collision at √ s NN = 5.02 TeV [108] is presented.

(3. 6 ) 3 . 1 .

 631 Fundamentals of Remler's Formalism 69 Proceeding with the time derivative inside the Tr[• • • ] and using the von Neumann equation (A.4), one gets

Chapter 3 .

 3 Remler Formalism: A Coalescence Model For Binding Particles based on the associated Wigner distribution W N .

FIGURE 3 . 1 -

 31 FIGURE 3.1 -Diagrammatic representation of Remler's effective rate by considering the change in the quarkonium formation probability during a collision

FIGURE 3 . 2 -

 32 FIGURE 3.2 -Wigner function width (Gaussian width) σ(T) as function of the temperature, obtained by solving the Schrödinger equation with a Rothkopf potential.The analysis of the dependence of the Gaussian width σ(T) with the temperature is restricted between the temperature values [0.15-0.4] GeV. Where T diss = 0.4 is the value of the dissociation temperature employed in our calculations.

.63) 3 . 5 .

 35 In-medium Q Q Interaction 87It can be used to eliminate the θ dependence of equation( 3.62)

FIGURE 3 . 3 -

 33 FIGURE 3.3 -Number of HQ pairs which relative distance in the C.M is lower than 1 fm as function of time for two different types of interaction as well as for free evolution.

Chapter 4 .

 4 Charmonium Production in pp Collisions thus, substituting in the equation (4.7) one gets

. 23 )

 23 Then from the set of equations (4.22) we set the components x r and y r of the vector r lab T as function of the r lab T,|| and r lab T,⊥ x lab r = r lab T,|| Cos(θ)r lab T,⊥ Sin(θ) y lab r = r lab T,|| Sin(θ) + r lab T,⊥ Cos(θ) . (4.24)

Chapter 4 .

 4 Charmonium Production in pp Collisionswill obtain the relation

  substituting the equation (4.27) into the (4.29) and integrating over the angle θ

Chapter 4 .

 4 Charmonium Production in pp Collisions probability one can obtain

FIGURE 4 . 1 -

 41 FIGURE 4.1 -Illustration of the dP J/Ψ dU T spectra behavior from the coalescence of c and c quarks under the approach A. In the result two different values of σ are presented together with three values of the initial separation (production separation) (r 0 ) as show in the legend of the figure

FIGURE 4 . 2 -

 42 FIGURE 4.2 -Illustration of the dP J/Ψ dU T spectra behavior from the coalescence of c and c quarks under the approach B. The same two different values of σ (0.303 and 0.35 fm) were consider, but now the initial separation (r 0 ) were consider in a distribution accordingly with the equation (4.43), in which four different values for the Gaussian width (σ r = 0, 0.044, 0.2 and 0.3 fm ) were consider

FIGURE 4 . 3 -

 43 FIGURE 4.3 -Comparison between the two approaches (A and B) for different values of the parameters, a direct value of r 0 in the case of approach A and r.m.s value (σ r ) of the distribution in the case of approach B.

(4. 52 )

 52 In the following two results, I present the theoretical prediction of our model for the ratio dσ J/ψ dy dσ c dy cfor the different values of interest of the free parameters (σ = 0.303 fm and σ = 0.35 fm) and r 0 (or σ r ) being r 0 = 0,0.044,0.2 and 0.3 fm. Each of the two following figures stand for each of the two approach (A and B) considered to take into account the initial separation.

FIGURE 4 . 4 -

 44 FIGURE 4.4 -Behavior of the

FIGURE 4 . 5 -

 45 FIGURE 4.5 -Same analysis, of the figure (4.4), but now considering σ = 0.35 fm

FIGURE 4 . 6 -

 46 FIGURE 4.6 -Comparison between our model J/Ψ production prediction for the different values of our model parameters (r 0 and σ) and the experimental data reported by ALICE collaboration in[119]

FIGURE 5 . 1 -

 51 FIGURE 5.1 -Comparison between our model prediction for the pp-scaled J/ψ direct production obtained from our Monte Carlo based code (blue line), together with the pp collision prediction obtained in the chapter 4 and the experimental results of ALICE collaboration [119] for the prompt production, where the scaled overlap function factor T PbPb ≈ 20.55 mb -1 represents an average value corresponding to the centrality bin from 0 -20%.

FIGURE 5 . 2 -

 52 FIGURE 5.2 -Comparison between our model prediction for the pp-scaled J/Ψ direct production obtained from our Monte Carlo based code (blue line), together with the experimental results from ALICE collaboration [124] scaled by the T PbPb ≈ 20.55 mb -1 . For the particular case of forward rapidity the number of c c pairs per unity of rapidity in the simulations are dN c c dy | y=0 = 13.12 and dN pp-scaled J/ψ dy

FIGURE 5 . 3 -

 53 FIGURE 5.3 -Time evolution of the J/ψ primordial production for various centrality range and considering two possible scenarios: one with the c c potential interaction ON (left) and other with the c c potential interaction OFF (right).

FIGURE 5 . 4 -

 54 FIGURE 5.4 -Time integral primordial p T distribution of J/ψ derived under four different evolution conditions combining c c potential interaction (ON and OFF) and medium elastic collision (ON and OFF).

FIGURE 5 . 5 -

 55 FIGURE 5.5 -Time evolution of the J/ψ primordial production obtained for interaction potential and medium elastic collision (ON) normalized by the square value of (( dN cc dy ) 2 ), for each centrality bin.

FIGURE 5 . 6 -

 56 FIGURE 5.6 -Number of c c pairs that at any given time fulfil the condition of r cm rel • p cm rel < h as function of total p T .

FIGURE 5 . 7 -

 57 FIGURE 5.7 -p T -spectra of the J/ψ integrated local rate production distribution for various centrality range and considering two possible scenarios: one with the c c potential interaction ON (left) and other with the c c potential interaction OFF (right).

FIGURE 5 . 8 -

 58 FIGURE 5.8 -Integral (time integral) local rate p T contribution to J/ψ production probability under four different evolution conditions c c potential interaction (ON and OFF) and medium elastic collision (ON and OFF).

FIGURE 5 . 9 -

 59 FIGURE 5.9 -p T -spectra of the average variation of the probability W J/φ coll due to a single interaction integrated over the full time evolution as indicated in the equation (5.8).

FIGURE 5 . 10 -

 510 FIGURE 5.10 -Integrated Collision rate p T -spectra obtained with c c interaction potential (ON) and medium elastic collision (ON) active and for three different collision criteria 0% (blue), 10% (red) as well as 5%(green) which is the usual value in the main calculations.

FIGURE 5 . 11 -

 511 FIGURE 5.11 -Distribution of momentum variation integrated over the whole evolution of all the c and c present in the evolution, one coming only from medium elastic collision (blue) and the momentum variation coming only for the potential (red).

FIGURE 5 .

 5 FIGURE 5.12 -Integrated Collision rate p T -spectra only due to c c interaction potential computed for three different Bjorken time steps (BTS): fullBTS (0.1 fm/c), halfBTS(0.05 fm/c) and quarterBTS(0.025 fm/c) and for our main Minkowski time step (0.25 fm/c).

FIGURE 5 . 13 -

 513 FIGURE 5.13 -Integrated Collision rate p T -spectra for a collision coefficient η cut = 5% for each of the possible combination of interaction potential and medium elastic ON and OFF, together with the trivial case OFF -OFF.

FIGURE 5 . 14 -

 514 FIGURE 5.14 -Time evolution of the J/ψ total production integrated over p T for different centrality bins and for two different scenarios: c c interaction potential ON (right plot) and OFF (left plot) and standard values for collision coefficient η cut = 5%, and Minkowski and Bjorken time steps 0.25 fm/c and 0.1 fm/c respectively.

FIGURE 5 .

 5 FIGURE 5.15 -Final J/ψ production probability p T -spectra obtained for different centrality bins and for two different scenario: Interaction potential OFF (left) and ON (right). The results were obtained with our usual parameter values for η cut = 5%, Minkowski time step value 0.25 fm/c and Bjorken time step 0.1 fm/c.

FIGURE 5 . 16 -

 516 FIGURE 5.16 -Final production probability p T -spectra for the four scenarios derived from interaction potential (ON and OFF) and medium elastic collision (ON and OFF), our standard values for collision coefficient η cut = 5%, and Minkowski and Bjorken time steps 0.25 fm/c and 0.1 fm/c respectively. In addition to our analysis, the ALICE experimental data[START_REF]Centrality and transverse momentum dependence of inclusive J/ψ production at mid-rapidity in Pb-Pb collisions at √ s = 5.02 TeV[END_REF] for mid-rapidity PbPb experimental data was including.

FIGURE 5 . 17 -

 517 FIGURE 5.17 -Transverse momentum distribution of the partial primordial production R AA (left) and partial integrated local rate R AA (right) of J/ψ production for the four different scenarios derived from medium elastic collision (ON and OFF) and interaction potential (ON and OFF).

FIGURE 5 .

 5 FIGURE 5.18 -Transverse momentum distribution of the partial integrated collision rate R AA (left) and the one corresponding to the final production (rigorous nuclear modification factor) R AA (right) of J/ψ production for the four different scenarios derived from medium elastic collision (ON and OFF) and interaction potential (ON and OFF).

FIGURE 5 . 19 -

 519 FIGURE 5.19 -Comparison between our model prediction for the nuclear modification factor R AA for two different centrality range 0 -20% (left) and 20 -40% (right) and the corresponding inclusive experimental data reported by ALICE [125] collaboration.

FIGURE 5 . 20 -

 520 FIGURE 5.20 -Nuclear modification factor behavior as function of the average number of participants in PbPb comparison versus the experimental data from ALICE (blue dots) [125], including the experimental error bars obtained by summing quadratically the statistical and systematic errors, obtained for interaction potential and medium elastic collision (ON) and our standard parameters values of: Bjorken and Minkwoski time step values of 0.1 fm/c and 0.25 fm/c respectively and η cut = 5%.

FIGURE 5 . 21 -

 521 FIGURE 5.21 -Comparison between our model prediction for the final "total" production (green) and also our model prediction for the diagonal contribution (red) for the case of medium elastic collision and interaction potential both ON, as well as our usual parameters value: collision coefficient η cut = 5% and Minkowski and Bjorken time steps 0.25 fm/c and 0.1 fm/c respectively.

FIGURE 5 .

 5 FIGURE 5.22 -Same comparison and with the same simulation parameters that the results presented in the figure A.3 but now for 20-40% (left) and 30-50%

FIGURE 5 . 24 -

 524 FIGURE 5.24 -Comparison between the equivalent nuclear modification factor of the diagonal contributions (red curve) and the off-diagonal contribution (blue curve) as a function of p T .

FIGURE 5 .

 5 FIGURE 5.25 -Equivalent Nuclear modification factor of the diagonal(left) and off-diagonal(right) component of the final production at 20-40% and 30-50% centrality bins as function of p T

FIGURE 5 . 26 -

 526 FIGURE 5.26 -Comparison between our model prediction for the mid-rapidity (|y| ≤ 0.9) elliptic flow v 2 for the centrality range 30 -50% and the experimental results reported by ALICE[START_REF]J/ψ elliptic and triangular flow in Pb-Pb collisions at √ s =5.02 TeV,ALICE Collaboration[END_REF] at mid-rapidity (left) and forward rapidity (2.5 < y < 4.0) (right), including the experimental error bars obtained by summing quadratically the statistical and systematic errors.

FIGURE 5 .

 5 FIGURE 5.27 -Comparison between our model prediction for the mid-rapidity (|y| ≤ 0.9) elliptic flow v 2 for the centrality range 10 -30% (left) together with 20 -40% (right) and the experimental results reported by ALICE[START_REF]J/ψ elliptic and triangular flow in Pb-Pb collisions at √ s =5.02 TeV,ALICE Collaboration[END_REF] at forward rapidity (2.5 < y < 4.0) including the experimental error bars obtained by summing quadratically the statistical and systematic errors.

FIGURE 5 .

 5 FIGURE 5.28 -Comparison between our model prediction for the mid-rapidity (|y| ≤ 0.9) elliptic flow v 2 , the elliptic flow of the diagonal component of the total production at final stage in the evolution and the experimental results reported by ALICE[START_REF]J/ψ elliptic and triangular flow in Pb-Pb collisions at √ s =5.02 TeV,ALICE Collaboration[END_REF] at forward rapidity (2.5 < y < 4.0) including the experimental error bars obtained by summing quadratically the statistical and systematic errors, for the centrality range 20 -40% (left) together with 30 -50% (right).

T 0 Γ

 0 Φ (t)dt avec Γ Φ défini comme Γ J/ψ (t) = dP Φ dt = d dt Tr[ρ Φ ρ N (t)]En appliquant la dérivée temporelle à l'intérieur de l'opérateur Tr [...] et en substituant l'équation de von Neumann, nous obtenonsP Φ (T) = Φ [H, ρ N ]]dt,En vue de simplifier les expressions, il s'avère intéressant de décomposer l'hamiltonien du système de N particules en trois termes:H N = H 1,2

W 2 σ 2 -p 2 σ 2 (hc) 2 où

 2222 Φ (r, p) = 8g (πh) 3 e -r g est un facteur de couleur de spin d'un méson vectoriel. Le facteur de spin du méson vecteur est3 4 et le facteur de couleur pour une paire Q Q de couleur neutre est vide (à 0 température), donc les effets liés à la température n'étaient pas considérés; 2-le potentiel d'interaction entre les particules qui forment l'état composé (particules 1 et 2) n'a pas été considéré, de sorte qu'au cours de leur évolution, les quarks lourds n'ont pas de moyen naturel d'empêcher leur éloignement progressif à l'intérieur d'un plasma en expansion; 3-aussi, pour que le modèle soit compétent aux énergies du LHC, il faut qu'il soit étendu aux domaines relativistes.

σ 2 (

 2

FIGURE A. 2 -

 2 FIGURE A.2 -Largeur de la fonction de Wigner (largeur gaussienne) σ(T) en fonction de la température, obtenue en résolvant l'équation de Schrödinger avec un potentiel de Rothkopf. L'analyse de la dépendance de la largeur gaussienne σ(T) avec la température est restreinte entre les valeurs de température [0,15-0,4] GeV où T diss = 0.4 GeV est la valeur de la température de dissociation utilisée dans nos calculs.

  Fig.A.3.

FIGURE A. 3 -

 3 FIGURE A.3 -Comparaison entre notre prédiction de modèle pour la production "totale" finale (vert) et aussi notre prédiction de modèle pourla contribution diagonale (rouge) pour le cas d'une collision élastique avec les particules du plasma et d'un potentiel d'interaction à la fois ON, ainsi que la valeurde nos paramètres habituels : coefficient de collision η cut = 5% et pas de temps de Minkowski et Bjorken 0.25 fm/c et 0.1 fm/c respectivement.

FIGURE A. 4 -

 4 FIGURE A.4 -Comparaison entre la prédiction de notre modèle pour le facteur de modification nucléaire de mid-rapidity (|y| ≤ 0, 9) R AA pour la plage de centralité 0 -20% (à gauche) avec la prédiction pour le flux elliptique v 2 pour la plage de centralité 20 -40% (à droite) et les résultats expérimentaux rapportés par ALICE [125] et [126] respectivement, y compris les barres d'erreur expérimentales obtenues par sommation quadratiquement les erreurs statistiques et systématiques.

  , T d = T n + T p and dN d dT d

	1.3. Coalescence Models of Quarkonium Production			27
	approximated as a function of the neutron and proton distributions as follows
	dN d dT d	=	p 3 0 6	m d m n m p	1 d + 2m d T d T 2	dN n dT n	dN p dT p	(1.16)
	where m d

TABLE 2 .

 2 1 -Summary of the lattice results on the dissociation temperature for Quarkonia states[START_REF]Heavy quarkonia in quark gluon plasma as open quantum systems[END_REF] and the vacuum binding energy values.

  is a 4-vector known as Minkowsky force and its spacial components are given by F i = K i /γ Chapter 3. Remler Formalism: A Coalescence Model For Binding Particles[115, 116]. The zero component one finds by the relation

  The different values P J/ψ correspond to the different combination of the parameter

	r 0	σ	0.303 fm	0.35 fm
	0 fm		0.101881	0.070423
	0.044 fm 0.09874 0.0687763
	0.2 fm	0.055698 0.0444154
	0.3 fm	0.028036 0.0260436

TABLE 4 .

 4 1 -Results of the conditional probability P φ obtained under the approach A, by integrating over U T the equation(4.45)

		Chapter 4. Charmonium Production in pp Collisions
	σ r	σ	0.303 fm	0.35 fm
	0 fm		0.101881	0.0705097
	0.044 fm 0.09888084 0.0688864
	0.2 fm	0.06366455 0.0483697
	0.3 fm	0.0439884 0.0352049

chosen in the previous results. Only significant differences between the approach A and B are observer for large values of r 0 ( or σ), meaning that those states with widest dimensions reflects the significant numerical difference.

TABLE 4 .

 4 2 -Results of the conditional probability P φ obtained under the approach B, by integrating over U T the equation(4.45)

TABLE 5 .

 5

1 -General values of average number of binary collisions (σ coll = 64 mb), the average number of c c pairs and total average number of primordial J/ψ in a typical collision in our simulations as function of the centrality bins.

  |p 1 |p 2 (A.[START_REF] Bodwin | Rigorous QCD analysis of inclusive annihilation and production of heavy quarkonium[END_REF] where we can identify Q with |κ| in the rest frame. Substituting the relative momentum into the expression of the state vector -Y))( f i (+|κ|)+ f i (-|κ|))|p 1 |p 2 -Y)) f i (+|κ|)|p 1 |p 2 (A.17)we have finally obtained the desired expression for |Ψ where the momentum p 1 and p 2 are taken as ( √ m 2 + κ 2 , ±κ) in the rest frame and then boosted with a rapidity shift +y B . Setting (√ m 2 + κ 2 , κ) = m(Cosh( ŷ), Sinh( ŷ), we get p 1 = m(Cosh(y B + ŷ), Sinh(y B + ŷ)) p 2 = m(Cosh(y Bŷ), Sinh(y Bŷ)).It remains to be checked that the orthogonality conditions of the state (A.17) is A.1. Finding Relativistic Orthogonal Two-Body Basis

	|Ψ y B ,i = fulfilled.	d|κ|dY √ s Φ y B ,i |Φ y B ,i = δ( √ sSinh(y B = dκdY √ s	δ(	√	sSinh(y B 163 dp 1 2E 1 dp 2 2E 2 δ( √ sSinh(y B -Y))δ( √ sSinh(y B -Y))×

  allows to create a completely orthogonal set of two-body particle states|Φ i = dy B g i (y B )|Φ y B ,i , (A.20)where g is an arbitrary function which satisfies the relationΦ i |Φ j = δ ij dy B g * i (y B )g j (y B ). (A.21)Realistically speaking, the best way to proceed is to work with the spectrumdN i dy B = | Φ y B ,i |Φ | 2 = |g i | 2 .

(A.22)

  + H N-2 + U 1,2 , où H 1,2 est l'hamiltonien des particules qui formeront l'état composé (Q et Q), H N-2 est l'hamiltonien du reste des N -2 particules du système et U 1,2 est le potentiel d'interaction entre les particules 1 et 2 avec le reste des particules du système, mais sans tenir compte de l'interaction entre 1 et 2.En appliquant la propriété cyclique de l'opérateur Tr [...] et en tenant compte du fait que [ρ Φ , H 1,2 ] = 0 car ρ Φ est l'opérateur densité issu des fonctions propres deH 1,2 et [ρ Φ , H N-2 ] = 0 car ρ Φ et H N-2 dépendent de différentes variables, nous obtenons : dP Φ (t) dt = Γ Φ (t) = -iTr[ρ Φ [U 1,2 , ρ N (t)]]Pour calculer le commutateur entre le potentiel d'interaction et l'opérateur de densité de probabilité, il faut se diriger vers l'espace des phases à l'aide des fonctions de Wigner correspondant respectivement à ρ Φ et ρ N Comme un traitement quantique complet de l'évolution de ρ N est hors de portée, on aura recours à une approximation semi-classique basée sur la distribution de Wigner associée W N . La densité de Wigner à N corps, W N , des quarks dans un plasma en expansion est un objet qui ne peut pas être calculé exactement. Dans le passé, il s'est avéré que de nombreuses observables dans les collisions d'ions lourds ont pu être bien décrites si l'on remplace la densité de Wigner à N corps par une moyenne sur les densités d'espace de phase classiquesW c W N ≈ W c N .Si les particules se dispersent seulement mais n'ont pas d'interaction potentielle, on parle d'approche en cascade ; si les particules interagissent en plus via un champ moyen ou par deux interactions corporelles on parle de "l'approche Boltzmann Uhling Uhlenbeck" et "Dynamique moléculaire quantique", respectivement. Dans les deux premières approches, La densité de Wigner classique est définie comme

	W c N = ∏

i h3 δ(r ir i (t))δ(p ip i (t))

  Le premier terme de l'expression Γ eff conduit au taux de collision, tandis que le second terme nous donne la variation de la probabilité de formation de quarkonium due aux variations de la densité de probabilité au cours du temps. En répétant le même processus utilisé pour obtenir le taux de collisions, c'est-à-dire en nous plaçant dans l'espace des phases via les fonctions de Wigner et en substituant ces fonctions dans le deuxième terme du taux effectif, nous obtenons La largeur gaussienne σ sera responsable de l'adaptation de la dépendance à la température du W Φ dans notre approximation. Cette dépendance sera dérivée de la relation entre σ et le rayon carré moyen de l'état lié (charmonium) et dont le comportement a la forme représentée sur la fig. A.2.

		r , T(t))	∂ρ N (t) ∂t	] + Tr[	∂ρ Φ (r, r , T(t)) ∂t	ρ N (t)]	(A.62)
	Γ loc = 8 Ṫ(t)	dσ(T) dT	∂ σ e	-( X 2 σ 2 + σ 2 p 2 h )
	=	16 (π) 3 σ(T(t))(	X 2 σ 3 (T)	-	σ(T)p 2 h2 )e	-( X 2 σ 2 + σ 2 p 2 h )

In practice, the potential is taken from recent work of[START_REF] Lafferty | Improved Gauss law model and in-medium heavy quarkonium at finite density and velocity[END_REF].

Mathematically, this is equivalent to take a vanishing Wigner density for T > T diss .

For consistency reason, the same potential[START_REF] Lafferty | Improved Gauss law model and in-medium heavy quarkonium at finite density and velocity[END_REF] was used for the Q Q semi-classical dynamics and the evaluation of the σ(T) width used in the Wigner distribution associated to the J/ψ(T) states.

Equal time t(Q) = t( Q) in the cm does not generally correspond to a unique time in the lab frame.

It is the overlap with the Wigner function that and the c c potential that generate the correlation.

Stemming from path-length asymmetry.

9 , qui conduit à g = 3 4 × 1 9 . Le facteur 8 (πh) 3 est liée à la normalisation de la densité de Wigner.En substituant les expressions des fonctions de Wigner W J/ψ et W c N dans l'équation

Chapter 5. Quarkonium production in AA collisions Nevertheless, comparing both scenarios (X1 and X2), we can notice that the ratio in X1 is less than in X2. In X1 the number of c c pairs per unity of rapidity is less (approximately half), so the number of off-diagonal pairs will be four times less than in the case X2, while the number of diagonal pairs will only be two times less than the case X2. Therefore, under the same hydrodynamic evolution, the diagonal contributions should increase their relevance with respect to the off-diagonal contributions the lower the number of c c pairs produced. Going deeper into our analysis, we compare now the results with the proton-proton collision scenario, which consists almost entirely of diagonal contributions. The enhancement at low p T of the experimentally observed production [START_REF]Centrality and transverse momentum dependence of inclusive J/ψ production at mid-rapidity in Pb-Pb collisions at √ s = 5.02 TeV[END_REF] in AA collision whit respect to the pp-scaled production can be explained from the point of view of the off-diagonal contributions.

Appendix A

General Formulation for Relativistic Wigner Function

The HQ are created in the early phase of the HIC through hard scattering processes and can have relativistic energies in the computational frame, the nucleus-nucleus center of mass frame. It is crucial that the Wigner function, used in our coalescence approach, is adequate to describe HQ under these conditions. We need therefore a relativistic invariant Wigner function which describes the dynamics that we are facing and remains boost invariant. It has also to give us the right expectation value for energy of a free particle

To obtain such a Wigner density we need first to build an orthogonal set of states that represents the bound states of Q Q and to ensure the correct normalization. From that we can then proceed to derive the Wigner function associated with the Q Q states, expressed in this orthogonal basis.

A.1 Finding Relativistic Orthogonal Two-Body Basis

Having this said, the first thing to do is to implement an orthogonal base on which we can project our two-body states. A quantum state for a two-body system (two scalar particles) can be written as

|Φ 2 and |Φ 1 are the Hilbert space vectors of the heavy quarks and |Ψ is the vector of the two body system. E 1 and E 2 are the total energies of the particles 1 and 2, respectively. p L,1 and p L,2 are the longitudinal components of the momenta and p T,1 and p T,2 the transverse momenta. The normalization of the single particles states is given by

Instead of normalizing the states like

we may consider final states with a well defined final center of mass momentum P, named Ψ P , normalized as

where P is the 3-momentum vector of the center of mass of the two-body system. We introduce states of a given center of mass momentum P by

where the index i refers to internal states. To present our approach we neglect the transverse degrees of freedom. This reduces our approach to one momentum space and one energy dimension. Then, substituting the later equation into the above normalization equation we obtain

We perform now a variable transformation to the center of mass and relative motion

] and get, assuming the center of mass motion is slow,

If the center of mass motion is slow the momenta of both particles have about the same absolute value and if the masses are identical we find

Finding Relativistic Orthogonal Two-Body Basis

161

. Requiring that f i and f j are orthogonal

one gets

which satisfies the energy normalization provided

with M = 2m Q . We boost now |Ψ from the center of mass system S to a system S' with p 1 = Λp 1 and p 2 = Λp 2 and consider the boosted state. Quite generically, one has

Because the measure (see eq. A.9) is invariant, we can write this as well as

To advance we introduce the invariant measure

and the relative momentum q = p 1 -p 2 2 = (q 0 , q) with d 2 p 1 d 2 p 2 = d 2 Pd 2 q. This yields

q 0 is given by the mass shell conditions δ(P • q) which yields

for |q| < P 0 /2, a condition which is systematically satisfied if the on-shell

Although we have now a method to build an orthogonal boost invariant base on which we can project our two-body states, the building of the Wigner function on this basis is not an obvious task. The main difficulty comes from the fact that it is not clear how to identify the conjugate variable of the rapidity. In the following sections of the appendix, we will discuss how a Wigner distribution can be obtained and for which states (states with well-defined momentum, rapidity, etc ... ) this is possible.

A.2 Implementation of a Relativistic Wigner Function from Given Basis in (1+1)D Case

The main difficulty for implementing a relativistic Wigner function from an orthogonal basis of the form (A.20), lies in the difficulty of defining a conjugate variable of the rapidity. We analyzed several conjugate variables for the rapidity. Among the options, only that to take the relative momentum p r in the center of mass as the conjugate variable allowed to arrive at a satisfactory conclusion and at the same time to obtain results which have clear physical significance. Here we discuss this option [Y, p r ].

We will show now how we obtain a Wigner density in [Y, p r ] for the basis of the form (A.20) To start, it is convenient to rewrite the (A.22) in the form

where the trace is performed over the phase space variables ([Y, p r ]), y Φ being the rapidity of the two body system (quarkonium). The density operators ρ Φ and ρ i,y Φ have the form Defining the relativistic velocity vector u = (u 0 , u) = (Cosh(y), Sinh(y)), we can write

we get

After introducing the auxiliary variables u = u+u 2 , p = p r +p r 2 , ∆u = uu and ∆p r = p rp r , we can write (A.29) formally as

We recall the relationship between conjugate variables

where x r is the relative position measured in the Q Q center of mass system and X is a dimensionless variable. We recall as well the definition of the Wigner function of a density operator ρ(r, r ):

where ρ(r, r ) = r|Φ Φ|r and thus substituting the delta functions relations (A. [START_REF]φ and Ω production in relativistic heavy-ion collisions in a dynamical quark coalescence model[END_REF] and comparing the definition of the Wigner function with the factors in the equation (A.30) we can derive the Wigner function for the two body density operator

and

Substituting the equations for the Wigner functions (A.34) and (A.33) in equation A. [START_REF] Gev | Production of Tritons, Deuterons, Nucleons, and Mesons by 30[END_REF] we obtain

Now we need to evaluate W i,y Φ (u, X; p r , x r ). For this we start from eq. A.17 and obtain:

.

. Substituting this expression in equation A. [START_REF] Greco | Quark coalescence for charmed mesons in ultrarelativistic heavy-ion collisions[END_REF] we obtain the desired Wigner density

With help of eq. A. [START_REF] Balk | Heavy Quark Effective Theory at Large Orders in 1/m[END_REF] we can write this expression as a function of the non relativistic wave function

2 )e -ix r ∆p r h is nothing than the non relativistic Wigner function W i,NR expressed in the coordinates of the center of mass. So we can finally establish

The equation for W i,y Φ does not depend on the variable X whose meaning was not very clear in physical terms (besides being the variable conjugate of ∆u). W i,y Φ depends on the rapidity of the center of mass of the quarkonia state in the computational frame. Equation (A.39) describes how we can evaluate the Wigner function for a Q Q pair in the bound state i with y Φ (u Φ ). It is also important to realize that the boosted Wigner function (A.39) inherits some properties of the non relativistic Wigner function W NR , including the normalization.

Substituting equation (A.39) in equation (A.30) we obtain the rapidity distribution of the Q Q pairs which are bound in a state i:

In equation (A.40) the Wigner function W i,NR (p r , x r ) represents the probability density (non-relativistic) of formation of a quarkonium state i. While W(u, X; p r , x r ) represents the probability density of finding a pair Q Q pairs with relative momentum and position q and x r and with center of mass velocity u. Since W references Q Q pairs before they form some bound state (free), the relationship with the two-body probability density operator ρ Q Q, can be trace back by using the Wigner function definition

where the integral over x r has the form :

which satisfies the following relationship:

is the total number of different Q Q pairs that can be formed. One can write the following relationship for W(u, p r , x r )

which is simply the differential form of equation (A.43). In our simulations, that magnitude will be interpreted as the quarkonia state density with a given 4-velocity and with a given relative momentum and relative position measured in quarkonia rest system. The rapidity distribution can be expressed as

A.3 Wigner Density in the Case of (1+1)D

It is necessary to evaluate the expression (A.46) of the previous Appendix in the coordinates of the quarkonia rest system for momenta given in the computational frame. We start to formulate the Wigner density for the case of (1 + 1)D; then it will be extended to the (3 + 1)D case. We define x CM and x lab as the relative coordinates measured in the quarkonium rest system and in the computational frame, respectively. We set x 0 lab = 0 because in the simulation programs we update positions and momenta with a fixed time step ∆t in the laboratory system. This implies that the relative time will be non-zero in the rest system of the Q Q pair. Under this condition the Lorentz transformation gives us:

where γ = u 0 . For the case (3 + 1) D this only applies to the longitudinal component:

In the center of mass system the energy of the Q and the Q are identical. If we define the 4-vector q = p 1 -p 2 2 we find as well that the zero component is zero in the center of mass system and we have

For the purposes of a Monte Carlo simulation, it is more natural to evaluate the magnitudes q y x in the quarkonium rest system. The rapidity distribution dN/dy of Quarkonia in a state I is therefore given by:

Appendix A. General Formulation for Relativistic Wigner Function

Note that we have:

and therefore we find

A.4 Generalizing for the (3+1)D Case

The results obtained in the case of (1 + 1) D, equation (A.46), can only be derived when one assumes a given velocity of center of mass (u). Therefore, it is only possible to derive an invariant Wigner function distribution function for those quarkonium states with a given velocity. The results derived in equations (A.45) and (A.46) can be generalized to the case (3 +1)D if we start from a state generated according to equation (A.2)

where (u 0 , u) is the 4-velocity of the quarkonium state Φ and p cm r is the relative momentum evaluated in the quarkonium state center of mass (cm). It is also possible to express these quarkonium states in the base (u T , Y), if we separate the velocity vector u into its transverse u T and longitudinal u z component, where Y = Atanh( u z u 0 ). Then generalizing equation (A.46) by extending this formula to the transverse component (3D case), we arrive at

A.4. Generalizing for the (3+1)D Case 171 where

(A.56) Equation (A.56) is the 3D generalization of equation (A.39). Inspection shows that we have just to multiply by the factor

. After performing the integrals over the delta functions we obtain:

This expression can be also expressed in the coordinates of the computational frame

in which the coordinates p cm r and x cm r have to be expressed as a function of the p lab r and x lab r . From the previous equation, integrating by the variables u T,Φ and y Φ we can obtain the absolute number of states of quarkonium i-th

In the Monte Carlo simulations, in which our model will be used, we take all possible (Q, Q) combinations into account in order to form i-th quarkonium state.

where the sum extends over all N Q × N Q possible combinations. So the total number of states Φ i is given by the expression

The normalization of the Wigner distribution function W i,NR should be normalized with the factor g = 2s Ψ +1 8(2s Q +1) 2 to take into account the spin of the quarkonium state and of the quarks. The factor 1/8 is due to the fact that only color neutral Q Q combinations can form a quarkonium cinétique, on obtient

où le facteur de la fonction delta de Dirac est égal à 1 si l'une des particules 1 ou 2 subit une collision avec une particule du milieu, comme le montre l'image A.1 . En revanche, s'il n'y a pas eu de collision pendant l'intervalle de temps donné, la fonction delta sera égale à 0.

Le deuxième facteur, composé de la projection du W J/ψ dans le W c N , peut être interprété comme la variation de la probabilité due à une collision avec les particules du milieu. Le taux Γ 1,2 indique la variation de la probabilité de formation d'un état composé (quarkonium) due aux collisions avec les particules du milieu, on l'appellera donc désormais taux de collision Γ 1,2 (t) = Γ coll (t).

Améliorations du modèle original de Remler

Le formalisme présenté ci-dessus doit être amélioré afin de décrire correctement la dynamique d'une collision d'ions lourds. Ces améliorations sont nécessaires car dans le formalisme original : 1-la formation d'états composés se produisait dans le Therefore, the construction of models that deal with the quarkonium dynamics and its time evolution is crucial to study the QGP matter. Despite that approaches like statistical hadronization rate equation models contributed greatly to the heavy-ion collision community. They have also pointed out that to get a completely satisfactory description of the problem of quarkonium production, a dynamical modeling of the recombination while preserving some quantum features and addressing the relativistic bound state formation, is need it.

To achieve that we develop a phase space coalescence model based on Remler's composite particle formation formalism to explore the quarkonium formation in a heavy-ion collision. Our model is a multiparticle-oriented formalism capable of dealing with the recombination probability of all possible heavy quarks pairs traveling in the QGP at a given time. By computing an effective differential rate, we can keep track of the probability through the whole evolution time as well as to study several observables that are probability depended such as the p T spectrum, R AA and flow coefficients (v 2 ). Additionally, these models offer an exciting window to address the inner dynamics of pair, including the interaction potential between the quarks and the influence of the medium on the pair dynamics.