.. |Tmunu| replace:: :math:`T^{\mu\nu}` .. |e| replace:: :math:`\epsilon` .. |u| replace:: :math:`u^\mu` .. |pimunu| replace:: :math:`\pi^{\mu\nu}` .. |Pi| replace:: :math:`\Pi` .. |Tmunu(X)| replace:: :math:`T^{\mu\nu}(X)` .. |e(X)| replace:: :math:`\epsilon(X)` .. |u(X)| replace:: :math:`u^\mu(X)` .. |pimunu(X)| replace:: :math:`\pi^{\mu\nu}(X)` .. |Pi(X)| replace:: :math:`\Pi(X)` .. role:: bash(code) :language: bash Before hydrodynamics: initial conditions, energy deposition & pre-equilibrium dynamics ======================================================================================= In these notes we define the early stage of heavy ion collisions as everything that happens before the plasma can be described by hydrodynamics. Hydrodynamics describes the *evolution* of the strongly-coupled quark-gluon plasma in space and time. It can describe this evolution if it is provided the *state* of the fluid on a given spacetime surface [*]_. .. [*] That is, if hydrodynamics is taken as nothing more than a system of partial differential equations, boundary conditions or initial conditions must be provided before a numerical solution of the equations can be obtained. .. There is no reason to believe hydrodynamics should describe the dynamics of heavy ion collisions at arbitrarily early time. Recall that, with second-order viscous hydrodynamics, the energy-momentum tensor |Tmunu(X)| is decomposed into four independent hydrodynamic fields: - The energy density |e(X)| - The flow velocity |u(X)| - The shear stress tensor |pimunu(X)| - The bulk pressure |Pi(X)| In almost every simulation of heavy ion collisions, the equations of hydrodynamics are initialized at a fixed value of :math:`\tau`\ : :math:`\tau_0`\ . That is, hydrodynamic initial conditions are provided as :math:`T^{\mu\nu}(\tau=\tau_0,x,y,\eta_s)`. More typically the different fields are initialized separately: :math:`\epsilon(\tau=\tau_0,x,y,\eta_s)`\ , :math:`u^\mu(\tau=\tau_0,x,y,\eta_s)`\ , :math:`\pi^{\mu\nu}(\tau=\tau_0,x,y,\eta_s)` and :math:`\Pi(\tau=\tau_0,x,y,\eta_s)`\ separately. Note that these initial conditions represent a considerable amount of information. To understand how much of a challenge it is to determine :math:`T^{\mu\nu}(\tau=\tau_0,x,y,\eta_s)`\ , it is worth taking a step back and look at the initial moments of a heavy ion collisions. The earliest stages of heavy ion collisions ------------------------------------------- No two collisions of nuclei are the same, and one of the reasons can be understood from a simple geometric picture, illustrated in the following figure: .. figure:: _static/impact_parameter_cern.jpg [Ref: http://cerncourier.com/cws/article/cern/53089 ] The two nuclei (ions) can overlap completely or can have a small overlap. The amount and the geometry of the plasma formed in these two cases will be very different. The geometric degree of overlap can be quantified with the *impact parameter* "b", shown on the above figure. In general, the more the overlap (the more "central" the collision is, "b" small), the larger the plasma. Smaller overlap between the nuclei results in smaller and shorter-lived plasma, and while very peripheral heavy ion collisions that barely graze each other might not produce any quark-gluon plasma at all. .. note :: The geometric picture discussed above is a very intuitive way of picturing heavy ion collisions. However, it is important to remember that the impact parameter "b" is not a quantity that can be extracted from heavy ion collision data. In this sense, it is preferable to think of it as a model parameter, not a measurable quantity. There are additional fluctuations aside from fluctuations in the degree of overlap of the nuclei: the plasma produced by two collisions with the same overlap will differ, since collisions of nuclei are quantum processes with a probabilistic rather than deterministic outcome. .. This implies that heavy ion collisions have a wide range of outcomes. .. .. topic:: Centrality classes Since heavy ion collisions can have different overlap and other fluctuations, they can produce A centrality class is a group of heavy ion collisions that have a similar multiplicity (the total number of hadrons measured by some detector). .. It is logical that collisions with similar impact parameters should be similar. It is thus very common to group .. use "proxys" for the impact parameter. For example, smaller impact parameters (more central collisions) should mean that a larger plasma is produced, which also implies that more final hadrons "\ :math:`N_\textrm{hadron}` \" are produced in such collisions. Thus, one could try to establish a relation between "\ :math:`N_\textrm{hadron}`\ " and "b". This can work well in general. The other major factor that affects the size and the lifetime of the plasma is the center-of-mass energy, :math:`\sqrt{s_{NN}}`. The higher the collision energy, the larger and longer-lived the produced quark-gluon plasma. .. note :: **MUSIC** and hydrodynamic models in general have been found to describe a large number of experimental measurements in collisions of gold nuclei (Au) at the top RHIC energy (:math:`\sqrt{s_{NN}}=200~\textrm{GeV}`) and collision of lead ions (Pb) at the LHC (:math:`\sqrt{s_{NN}}=2760~\textrm{, } 5020~\textrm{GeV}`, ...). It is likely that hydrodynamic models can still be used in collisions of smaller nuclei and/or lower :math:`\sqrt{s_{NN}}` collisions, but one should proceed with additional caution as these systems can have additional complications. Obtaining initial conditions for hydrodynamics ------------------------------------------------ From a purely practical point of view, we can ask the question "Given the species of the colliding nuclei (e.g. Au or Pb) and a collision energy (\ :math:`\sqrt{s_{NN}}`\ ), what are the initial conditions :math:`T^{\mu\nu}(\tau=\tau_0,x,y,\eta_s)` of hydrodynamics?" This is a pragmatic question that hides a lot of the complexity of the early stage of heavy ion collisions. Between the moment that the nuclei collide and the time at which hydrodynamics become applicable, there is a complex dynamical evolution of the deconfined matter. .. Considerable efforts are being made to better understand the transition from the Moreover the fact that there are fluctuations in heavy ion collisions implies that there is no single initial condition corresponding to two nuclei colliding. There is rather a family, or a distribution, of initial conditions. Describing this pre-hydrodynamic stage of stage of heavy ion is out of the scope of **MUSIC**. Older versions of MUSIC included simple models of initial conditions such as Optical Glauber and Monte Carlo Glauber. This is not the case anymore: MUSIC must now be provided with an external initial condition, such as Trento, IP-Glasma, etc. .. Broadly speaking, there are two options to initialise the hydrodynamics in **MUSIC** at :math:`\tau=\tau_0`: .. .. - An external model can be used to describe the early stage of to collision. The results of this external model can be read in and used as initial conditions for **MUSIC**. This is the procedure adopted when IP-Glasma initial conditions are used with **MUSIC**. .. .. - The initial conditions of hydrodynamics can be parametrized with the Glauber ansatz .. Considerable progress has been made to describe this transition to equilibrium from first principles, but in general it is still an unsolved problem. The initialization of the hydrodynamics fields is performed in function "Init::InitTJb()" in file `"init.cpp" `_. This function calls various subroutine meant to read external initial condition files in different formats. .. Please refer to these subroutines for additional information on the format. .. The input parameter :bash:`Initial_profile` is used to choose which type of initial conditions are used to the subsequent hydrodynamic evolution. The most important options are .. .. code:: .. .. ################################### .. # parameters for initial conditions .. ################################### .. Initial_profile 3 # type of initial condition .. # 1: Optical Glauber model .. # 3: MC-Glauber model .. # 8: Read in initial profile from a file .. # (e.g., IP-Glasma) .. .. .. Reading-in an external initial condition is a more flexible option, but in certain cases Glauber initial conditions are sufficient. Reading-in external hydrodynamic initial conditions ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Since there is no accepted universal format of hydrodynamic initial conditions, reading-in external initial conditions will almost certainly require modifications to the code. This should be straightforward, however, given examples already in **MUSIC** that show how to initialize hydrodynamics with external files. The name of the external file containing the initial conditions is specified by the input parameter, :bash:`Initial_Distribution_Filename`. The expect format for this initial condition file is communicated with the input parameter :bash:`Initial_profile` For example, the file format expected when using :bash:`Initial_profile 8` is defined in the subroutine "Init::initial_IPGlasma_XY()". First, there is a header which indicates the number of grid points in each direction and the size of the grid cells, in fermi for x/y (unitless in :math:`\eta`\ ): .. code-block:: C++ string dummy; double dummy2; profile >> dummy >> dummy >> dummy2 >> dummy >> neta >> dummy >> nx >> dummy >> ny >> dummy >> deta >> dummy >> dx >> dummy >> dy; where "neta", "nx" and "ny" are the number of grid points in each directions, and "deta", "dx" and "dy" are the step size. An example of header line would be: .. code:: # dummy 1 etamax= 1 xmax= 200 ymax= 200 deta= 0 dx= 0.17 dy= 0.17 The rest of the file is read as follows: .. code-block :: C++ for (ix = 0; ix <= DATA->nx; ix++) { for (iy = 0; iy <= DATA->ny; iy++) { profile >> dummy1 >> dummy2 >> dummy3 >> density >> utau >> ux >> uy >> dummy >> dummy >> dummy >> dummy; where "density" is the energy density in GeV/fm^3, utau is u^\tau, ux is u^x and uy is u^y. The "dummy" variable contain either redundant or unused information. The information read from the external file is then used to initialize the discretized hydrodynamic fields: .. code-block :: c++ void Init::initial_IPGlasma_XY(int ieta, SCGrid &arena_prev, SCGrid &arena_current) { [...] // read the external initial condition file double density, dummy1, dummy2, dummy3; double ux, uy, utau; for (int ix = 0; ix < nx; ix++) { for (int iy = 0; iy < ny; iy++) { profile >> dummy1 >> dummy2 >> dummy3 >> density >> utau >> ux >> uy >> dummy >> dummy >> dummy >> dummy; temp_profile_ed[ix][iy] = density; temp_profile_ux[ix][iy] = ux; temp_profile_uy[ix][iy] = uy; temp_profile_utau[ix][iy] = sqrt(1. + ux*ux + uy*uy); [...] } } profile.close(); // Initialize the different component of T^{\mu\nu} with this information double eta = (DATA.delta_eta)*ieta - (DATA.eta_size)/2.0; double eta_envelop_ed = eta_profile_normalisation(eta); int entropy_flag = DATA.initializeEntropy; for (int ix = 0; ix < nx; ix++) { for (int iy = 0; iy< ny; iy++) { double rhob = 0.0; double epsilon = 0.0; if (entropy_flag == 0) { epsilon = (temp_profile_ed[ix][iy]*eta_envelop_ed *DATA.sFactor/hbarc); // 1/fm^4 } else { double local_sd = (temp_profile_ed[ix][iy]*DATA.sFactor *eta_envelop_ed); epsilon = eos.get_s2e(local_sd, rhob); } if (epsilon < 0.00000000001) epsilon = 0.00000000001; arena_current(ix, iy, ieta).epsilon = epsilon; arena_current(ix, iy, ieta).rhob = rhob; arena_current(ix, iy, ieta).u[0] = temp_profile_utau[ix][iy]; arena_current(ix, iy, ieta).u[1] = temp_profile_ux[ix][iy]; arena_current(ix, iy, ieta).u[2] = temp_profile_uy[ix][iy]; arena_current(ix, iy, ieta).u[3] = 0.0; arena_prev(ix, iy, ieta) = arena_current(ix, iy, ieta); } } } To summarize, there are two main input parameters used in this case: .. code :: bash Initial_profile 8 # 8: Read in initial profile from a file # (e.g., IP-Glasma) # external file from which the initial conditions are read Initial_Distribution_input_filename example_inputs/filename .. topic :: Extra parameters For the example above, there are two extra parameters. - Parameter :bash:`initialize_with_entropy` should be used if the initial condition file provides entropy density instead of energy density. This is an unlikely occurence and should likely not be used. - Parameter :bash:`s_factor` is designed to be used to normalize the energy density read-in from the external file :bash:`Initial_Distribution_input_filename`. Most of the time this parameter should not be used and should be set to "1.0". .. code-block :: bash initialize_with_entropy 0 # 0: with energy density # 1: with entropy density s_factor 1.0 # normalization factor read in # initial data file .. Parameter :bash:`s_factor` is designed to be used to normalize the energy density read-in from the external file :bash:`Initial_Distribution_input_filename`. Glauber initial conditions ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Older versions of MUSIC included Glauber initial conditions. This is no longer the case. .. .. A comprehensive presentation of the Glauber model can be found in e.g. https://inspirehep.net/record/742696\ . There are two main inputs to the Glauber model: the distribution of nucleons inside the nuclei, and the inelastic nucleon-nucleon cross-section. In very simplified terms, if the nuclei colliding in heavy ion collisions are taken as bags of nucleons, the Glauber model describes the collision of two nuclei as the incoherent collision of the all the nucleons inside the nuclei. This is illustrated in the following figure. .. .. .. figure:: _static/mc_glauber_holopainen.png .. .. [Ref: `Hannu Holopainen, PhD thesis `_ ] .. .. Every dark circle represents a nucleon that underwent a nucleon-nucleon collisions. They are referred to as "participant nucleons". Light circles are nucleons that did not collide with any nucleons from the other nuclei. Nucleon-nucleon collisions are called binary collisions. .. .. The Glauber model can be used as initial condition for the energy density |e(X)| in the plane transverse to the collision axis (x-y plane). Energy can be deposited either at the position of the binary nucleon collisions, or at the position of the participant nucleons, or a combination of the two. An input parameter in **MUSIC** controls :bash:`binary_collision_scaling_fraction` the fraction of energy density initialized with binary collision and participants nucleons. .. .. The picture above is technically a *Monte Carlo* implementation of Glauber, referred to as Monte Carlo Glauber or MC Glauber. There is another version of Glauber, *optical* Glauber, which is roughly equivalent to averaging a large number of MC Glauber initial conditions, producing a smooth energy density distribution. Both MC Glauber and optical Glauber are available in **MUSIC**. .. .. Here is a summary of which parameters are used for optical Glauber, which is used for MC Glauber, and which is used for both: .. .. .. code :: .. .. ################################### .. # parameters for initial conditions .. ################################### .. Initial_profile 3 # type of initial condition .. # 0: for Gubser flow test (future) .. # 1: Optical Glauber model .. # 3: MC-Glauber model .. # 8: Read in initial profile from a file .. # (e.g., IP-Glasma) .. # .. initialize_with_entropy 0 # 0: with energy density .. # 1: with entropy density .. # .. s_factor 30.0 # normalization factor for the energy/entropy .. # applied to MC-Glauber initial conditions .. # (Initial_profile 3) .. # or initial data file (Initial_profile 8) .. # .. Eta_plateau_size 20. # size of the plateau in spatial rapidity eta_s .. Eta_fall_off 0.7 # the scale of the fall off of the .. # plateau in eta_s direction .. .. # .. # parameters for Glauber model (Optical and MC Glauber, Initial_profile == 1 or 3) .. binary_collision_scaling_fraction 0.14 # for wounded nucleon/binary .. # collision mixing ratio in the .. # Glauber model .. SigmaNN 42.1 # NN inelastic cross section for .. # Glauber model .. Target Au # type of target nucleus .. Projectile Au # type of projectile nucleus .. .. # parameters for optical Glauber model only (effective only when Initial_profile == 1) .. Impact_parameter 9. # impact parameter (fm) .. # .. Maximum_energy_density 54. # maximum energy density for Glauber .. # initial condition .. # (used instead of parameter "s_factor") .. .. # parameters for MC Glauber model only (effective only when Initial_profile == 3) .. sigma_0 0.4 # the width of the gaussian .. # for entropy density profile .. bmin 0.0 # minimum value of impact parameter .. bmax 2.0 # maximum value of impact parameter .. .. .. .. topic :: Note on the Glauber implementation in **MUSIC** .. .. The implementation of MC Glauber used in MUSIC is described in Section IV of https://arxiv.org/pdf/1109.6289.pdf\ . This implementation may differ from other implementations of MC-Glauber, in particular in the following aspects: .. .. - There is no minimum distance between nucleon sampled from the nucleus .. - Negative binomial fluctuations are not included in the code, that is all nucleon-nucleon collisions deposit the same amount of energy/entropy .. - Energy/entropy is deposited at the center of mass of binary collisions. .. .. .. .. Parameters "Target", "Projectile", "SigmaNN", "bmin" and "bmax" should map without any ambiguity to any other implementation of MC-Glauber. .. .. .. Parameter "binary_collision_scaling_fraction" correspond to a fairly standard procedure in the field, there shouldn't be any ambiguity about it either. .. .. .. Parameter "sigma_0" is the spatial width of the Gaussian used to deposit energy/entropy, as in :math:`\exp\left(-r^2/2/(\textrm{sigma_0})^2\right)` There might be slight differences in the definition of "sigma_0" in other Glauber implementation. .. .. Since the Glauber model only provides an initialization for the energy density |e(X)| in the transverse plane, it must be complemented by other models or assumptions to provide a complete initial conditions. To obtain a profile in the :math:`\eta_s`\ , the input parameters :bash:`Eta_plateau_size` and :bash:`Eta_fall_off` can be used to multiply the transverse energy profile by a function that has a plateau around midrapidity (\ :math:`\eta_s \lesssim` :bash:`Eta_plateau_size` ) and falls off with a width :bash:`Eta_fall_off` after. .. .. This still only provides an initialization for the energy density |e(X)|. When Glauber initial conditions are used, **MUSIC** initialises the flow velocity |u(x)| to :math:`(1,0,0,0)` and the viscous parts of the energy-momentum tensor (|pimunu(X)| and |Pi(X)|) to zero. This can only be modified by modifying the code in "init.cpp".