Holographic heavy ion collisions with baryon charge

We numerically simulate collisions of charged shockwaves in Einstein-Maxwell theory in anti-de Sitter space as a toy model of heavy ion collisions with non-zero baryon charge. The stress tensor and the baryon current become well described by charged hydrodynamics at roughly the same time. The effect of the charge density on generic observables is typically no larger than 15\%. %The rapidity profile of the charge is wider than the profile of the local energy density. We find significant stopping of the baryon charge and compare our results with those in heavy ion collision experiments.


Introduction
Holography has provided useful toy models with which to study the quark-gluon plasma (QGP) created in heavy ion collision experiments. In these models one considers a caricature of the incoming ions consisting of two gravitational shock waves moving at the speed of light towards each other in anti-de Sitter space (AdS) [1][2][3][4]. The collision results in the formation and subsequent hydrodynamization of a strongly coupled QGP which can be studied from first principles (see e.g. [5,6] for recent reviews).
The dynamics of the gravitational field in AdS encodes only the dynamics of the stress tensor in dual gauge theory. In other words, the collisions considered in the references above result in the formation of a neutral plasma. However, the plasma created in heavy ion collisions certainly carries some non-zero baryon charge which increases in importance as the energy of the collision decreases. The goal of this paper is therefore to provide the first simulation of charged shock wave collisions in AdS in order to model the formation and hydrodynamization of a QGP with a non-zero baryon charge. The conserved baryon current is dual on the gravity side to a Maxwell field, so we consider collisions in Einstein-Maxwell theory. 1

The model
We start with the Einstein-Maxwell action with a negative cosmological constant where G is Newton's constant, R is the Ricci scalar, L is the asymptotic AdS radius, e is a parameter controlling the backreaction of the Maxwell field on the metric, and F mn ≡ ∂ [m A n] is the electromagnetic field strength with A m the vector potential. The metric and the gauge field are respectively dual on the gauge theory side to the stress tensor T µν and to a conserved U (1) current J µ whose time component one is free to think of as the baryon number density. In the case in which the action (2.1) is viewed as a consistent truncation of the dimensional reduction of type IIB supergravity on S 5 the dual gauge theory is N = 4 Super Yang-Mills and the U (1) current arises from the R-symmetry of this theory. 2 When we need to be concrete (for example to fix normalization factors) we will adopt this viewpoint. However, we emphasize that for most purposes the specific origin of the Maxwell field is unimportant and one could think of (2.1) simply as a bottom-up model that incorporates the minimal set of ingredients to describe the dynamics of the stress tensor and a conserved U (1) current in the dual gauge theory.
The equations of motion following from (2.1) are is the stress tensor sourced by the electromagnetic field. These equations admit the following solution describing a charged shock wave moving at the speed of light: where x ± = t ± z and h(x + ) and a(x + ) are arbitrary functions of x + . In the case in which the dual theory is N = 4 SYM these functions are related to the expectation values of the corresponding dual operators through [9] T ++ = N 2 The fact that J scales as N 2 c reflects its R-symmetry origin, namely that microscopically it is built out of adjoint degrees of freedom. An analogous solution describing a wave moving in the opposite direction is obtained by replacing x + → x − in the expressions above.
We will adopt the following choices for the functions h and a describing our incoming projectiles: The energy and charge densities per unit area of the shock are 2π 2 m 3 /N 2 c and 2π 2 m 2 /N 2 c , respectively. Note that we choose the centers and the widths of the Gaussian profiles for h and a to be the same, as corresponds to the fact that we want to model the collision of projectiles that carry both energy and charge.
We see from Eq. (2.2) that the parameter e controls the magnitude of the backreaction of the Maxwell field on the dynamics of the spacetime metric. It is clear from the action (2.1) that this parameter could be absorbed in the normalization of the Maxwell field, at the expense of including an explicit factor of e multiplying the gauge field amplitude a(x + ) in Eq. (2.6). In other words, one can either think of the backreaction as being controlled by e for a fixed incoming amplitude a(x + ), or as it being controlled by the incoming amplitude for a fixed e. We find it convenient to adopt the first viewpoint.
The limit e → 0 on the gravity side corresponds to an approximation in which the Maxwell field is treated as a probe field that propagates in a fixed background without affecting it. In the gauge theory this means that the charge density is treated in a quenched approximation in which it has a negligible effect on the dynamics of the gluons. This limit is physically interesting since it describes the situation in which the energy density dominates over the charge density, and also as a benchmark against which the backreacted results can be compared. In Sec. 4 we will consider the collision in the probe approximation, and in Sec. 5 we will consider the backreacted case.

Thermodynamics and hydrodynamics
We begin by reviewing the thermodynamics of the plasma in the theory described by the action (2.1). This will be useful when we later consider the hydrodynamics of the plasma.
We will work with a rescaled stress tensor and current given by and set The relation between these quantities and the temperature T and the chemical potential µ is given by the following expressions (see e.g. [10]): Using these expressions one can show that, for a small charge density, the chemical potential is given by (3.6) We see that, for fixed energy density, µ and ρ are linear in one another with a 1% (10%) accuracy if ρ/E 3/4 is no larger than 0.15 (0.48). In the opposite limit, when the chemical potential is high compared to the temperature, y 1, we have that and the ratio of charge to energy density approaches the value The subindex 'max' indicates that this is the maximum value of ρ for a fixed energy density. At this value the black brane becomes extremal and the temperature approaches zero compared to any other scale. Above this value a naked singularity appears. These features are illustrated in Fig. 1.
Charged hydrodynamics for N = 4 SYM has been studied extensively [10,11]. In the case of a 1+1 dimensional flow, as is of interest here, the constitutive relations take the form 3 where (3.12) and the transport coefficients (rescaled by (3.14)

The probe approximation
With the initial conditions (2.9) for the collision we can solve the Einstein-Maxwell equations, as outlined in [12,13]. 4 In this section we work in the strict probe approximation for the Maxwell field. We are therefore computing the evolution of a Maxwell field on top of the dynamical shock wave collisions studied earlier in [1][2][3][12][13][14]. Computationally it is however convenient to evolve both the metric and the Maxwell field at the same time, thereby recomputing the gravitational shock wave background.
We performed simulations for mw = 0.1 (1/4 -shocks in the language of [2]) and mw = 1.9 (2 -shocks in the language of [2]), which we will refer to as thin and thick shock collisions, respectively. The resulting charge density ρ is plotted in Fig. 2 for both cases, whereby we included the energy density E for comparison (also found in [2]). Clearly in the thin regime the shocks gradually lose their charge into a charged plasma between the shocks, much like the energy density. For the thick shocks the energy density already hydrodynamizes during the collision regime [2] and, as we will see below, the same is true for the charge density. The shape of the charge and energy densities is different, a feature which we will later analyze in the local rest frame as well.
One interesting feature is the decay of the original shocks in the thin regime. Indeed, we see that on the attenuating maxima the charge density exactly follows the energy density, despite the fact that the charge and the energy are distributed differently in between the shocks. This agreement on the attenuating maxima is remarkable in view of the fact that the energy and the charge are governed by in principle completely different dynamics, i.e. the Einstein equations and the Maxwell equations on a fixed background, respectively. This suggests that on the light cone a simplified picture may be possible (perhaps along the lines of [15]). Having the complete stress-tensor we are able to extract the energy and charge density in the local rest frame, defined through from which we find where u µ is the (timelike) fluid velocity. For the charge density we have This allows us to study the approach to charged hydrodynamics. We read off from our simulations the fluid velocity and the energy and charge densities. We then use the constitutive relations (3.9)-(3.10) to obtain the hydrodynamic prediction for the transverse and longitudinal pressures, P hyd L and P hyd T , and for the time component of the current, J hyd t . We define the hydrodynamization time for the stress tensor, t hyd , as the time beyond which Similarly, we define the hydrodynamization for the current, t J hyd , as the time beyond which A crucial difference between the hydrodynamization of the stress tensor versus the charge current is that the pressures at z = 0 can deviate from their hydrodynamic values due to the gradients of the velocity field. In contrast, parity symmetry implies that u has no spatial component at z = 0 and that J z vanishes identically. As a consequence, the (ideal) hydrodynamic prediction for the current is always exact at mid-rapidity. In order to asses the validity of the hydrodynamic description for J we will therefore look at non-zero rapidity.
The hydrodynamization times for the stress tensor in the case of a neutral fluid were determined in Refs. [1,2,14]. The result is 5 Since in this section we are working in the probe approximation, these times agree with our hydrodynamization times. With our new simulations we can now also study the hydrodynamization of the charge current, whose evolution and hydrodynamic approximations are shown in Fig. 3. The hydrodynamization times are In the case of thick shocks we have listed the value t J hyd = 0 to indicate that the current is always well predicted by hydrodynamics. In contrast, for thin shocks we see that the hydrodynamization times for the current away from mid-rapidity are very similar to those for the stress tensor. In both cases, after t J hyd there is a small but significant difference between ideal and viscous hydrodynamics, which shows that the fluid velocity of the charge differs from the velocity of the energy-momentum flow. Fig. 4 shows the spacetime rapidity profile of the local charge density. This is compared with the local energy density rapidity profile for several different proper times in Fig. 5.
Especially striking in Fig. 5 is the development of maxima at non-zero rapidity for the charge profile; we will come back to this in the Discussion section. This happens fast in the Figure 3. Time component of the current J t (black, solid curve) and its approximations based on ideal (blue, solid curve) and viscous (blue, dashed curve) hydrodynamics, for thin (top row) and thick (bottom row) shocks. The hydro curves start at a time after which a local rest frame can be defined (as in [16]). The insets show that even though ideal hydrodynamics gives a better overall fit in the range plotted, viscous hydrodynamics gives a better description of the final approach to hydrodynamics. Thick shocks (bottom) are always well described by charged hydrodynamics, whereby viscous hydrodynamics gives a significant improvement.  [16], and the hydrodynamic region is indicated by the dashed line (according to (4.4)). The upper white regions are outside our numerical grid. The right plot is for thick shocks, which is always in the hydrodynamic regime for τ > 2. The black lines are stream lines of the fluid velocity, which are similar to the charge velocity. Figure 5. Rapidity profile of the local charge (solid curves) and the local energy (dashed curves) densities at different proper times τ as a function of rapidity y and normalized at mid-rapidity (see also [2,13,14] for similar plots of the energy density). For thin shocks the energy density does not evolve much on the times shown, whereas the profile for the charge density widens significantly. For thick shocks the charge and energy profiles widen in a similar fashion, though the profile of the charge is always wider. case of thin shocks, whereas it takes longer for thick shocks. Furthermore, we see that for thin shocks the evolution of the charged rapidity profile is expanding much faster than the equivalent profile for the local energy density. For thick shocks both profiles expand at a similar rate.
To characterize the charge deposition after the collision, we compare the charge density in the local rest frame, ρ loc , to the local energy density, E loc . More precisely, at all points in spacetime at which a local rest frame exists, we define the ratio The significance of this ratio, which we plot in Fig. 6, is that, in equilibrium, µ eff coincides with the small-charge limit of chemical potential of the plasma, Eqn. (3.6). We will come back to Fig. 6 in the Discussion section. When colliding shocks with symmetric energy and charge distributions it is impossible to determine which part of the energy and charge in the final plasma comes from the leftor right-moving shock. In our current set-up, however, it is possible to charge only the leftmoving shock, leaving the right-moving shock neutral. This results in a charge distribution as shown in Fig. 7. 6 Clearly most of the charge ends up at z < 0, but at later times a surprisingly large fraction of about 41% ends up at z > 0, indicating that the strong interactions of the collision can indeed let the charge density bounce back, reflecting the direction of 41% of the charge.
Quite interestingly, the right-moving charge contains a bump moving close to the speed of light, though not as fast as the left-moving shock moving at the speed of light. This bump, however, does not come with a minimum with negative charge density, indicating that the negative charge density has to be associated to the presence of the original shock on the light cone. We verified that the profile away from the original shock is independent of the width of the charge, in agreement with the findings in [3]. Experimentally, these simulations have potential consequences when it is possible to vary the baryon charge while leaving the energy constant. This can for instance be done by comparing proton with antiproton collisions, or deuteron collisions with proton collisions of double the energy (preferably selected such that the deuteron is aligned along the beam direction).

Charge backreaction
We now move away from the case e = 0 and include the full backreaction of the charge density. We find that, for a given width of the incoming shocks, there is a maximum value of e that our code is able to evolve. One clear physical reason why there is such a maximum value is that the charged plasma formed after the collision has to have a smaller charge density than the maximum density set by Eq. (3.8), and indeed our simulations reach as high as 80% of the maximum value according to the instantaneous energy density. On the other hand, this maximum value is also sensitive to the performance of the numerical code and furthermore depends sensitively on the initial separation of the shocks. The latter may be an indication that single charged shock waves cannot be evolved in a stable manner by themselves. Nevertheless, choosing initial conditions with the shocks close enough to each other determines well-defined initial conditions, which shows us the effect of the back reaction of the charge on the metric. In a nutshell, the conclusions from these simulations are that (i) the effect of the charge on generic observables is relatively small, and (ii) this effect scales approximately linearly with e 2 . We emphasize that these results are non-trivial since the charge density attains values well into what is expected to be the non-linear regime based on Eq. (3.6).
To illustrate these results, in Fig. 8 we compare the energy density at mid-rapidity as a function of time for collisions with identical initial conditions but different values of e. The first row shows the energy density itself, whereas the second row shows the relative difference between a given collision and the e = 0 collision, normalized by e 2 . The maximum value of the curves on the second row tells us that the maximum effect of the backreaction is about 0.12 × 1.12 2 15% for thin shocks and 0.0075 × 2.85 2 6% for thick shocks. The fact that the curves on the second row fall almost on top of each other for thin shocks (left) means that in this case the backreaction of the charge is almost exactly linear in e 2 . For thick shocks the deviations from linearity are slightly larger. Note that in both cases the charge density as a function of time is not small. Indeed, for thin shocks it initially exceeds 60% of the would-be maximum value according to the instantaneous energy density and drops to around 30% at later times. For thick shocks the maximum exceeds 80% and the curve drops to 70% at later times. Fig. 9 shows the same information as Fig. 8 but as a function of z at a constant time t 0 . In the case of thin shocks we have chosen t 0 = t hyd , whereas in the case of thick shocks t 0 is the time at which the green curve in the second row of Fig. 8 attains its maximum.
As expected from this discussion, the effect of the charge density on the hydrodynamization time of the stress tensor is relatively small. For example, we find that for collisions with mw = 0.1 the hydrodynamization time in the charged case, with e = 1.07, is 3% shorter than in the neutral case. Instead, for collisions with mw = 0.75 the hydrodynamization time in the charged case, with e = 1.7, is 6% longer than in the neutral case.

Discussion
We have studied collisions of charged shocks in AdS 5 . Via the gauge/gravity duality, these are dual in the gauge theory to collisions of lumps of energy carrying fixed amounts of baryon charge (per unit area). As discussed in [2], the dynamics of the stress tensor show qualitatively different features depending on the width of the colliding shocks. In this paper we have shown that similar qualitative differences also appear in the distribution of baryon charge after the collision. One of the main observations of [2] was that, while thick shocks lead to a complete stopping of the incident energy followed by a subsequent hydrodynamic evolution, narrow shocks exhibit a transparent regime at early times in which the initial shocks cross each other depositing their energy gradually as time progresses. In this paper we have shown that the baryon charge deposition exhibits identical behaviour, as illustrated in Fig. 2.
This transparent regime is transient, meaning that at a sufficiently late time the receding shock fragments are completely absorbed by the plasma. Nevertheless, the space-time rapidity distribution of baryon charge exhibits interesting diverse features for the two cases we have explored, as shown in Fig. 4. At fixed proper times, both for thick and thin shocks the charge distribution is wider in rapidity than the energy distribution. This means that, as rapidity grows, the plasma becomes more baryon rich, as illustrated by the increase in the chemical potential plotted in Fig. 6. However, for thin shocks the deposition of charge is wider than for thick shocks and it evolves much faster with time. Quite remarkably, for thin shocks the initial close-to-Gaussian charge distribution at early times, prior to hydrodynamisation, changes quickly to an almost-flat distribution in space-time rapidity which, within our limited numeric range, hints at the formation of maxima at relatively large rapidity, y > 1. In the range of proper times and rapidities covered by our simulations, see Fig. 5, these structures appear in regions in which the evolution is well described by hydrodynamics. However, the formation of these maxima may involve far-from-equilibrium dynamics that the hydrodynamic regions are in causal contact with. For thin shocks at very late times, the large rapidity region also develops local minima close to the edge of the rapidity coverage, which arise solely from the hydrodynamic evolution of the plasma.
The space-time rapidity profile of the charge distribution hence indicates that collisions of shock waves lead to a significant stopping of the baryon charge of the incident projectiles. To best illustrate this point it is instructive to determine the fraction of the total charge (per unit area) of the incident shocks that is deposited between momentum rapidity y p = −1 and y p = 1, with y p = arctanh(v z ) and v z the velocity field in the collision direction. This is determined by integrating the charge between the two points where the fluid velocity Figure 10. We show the fraction of the charge (given in Eq. (3.2)) in the plasma that has momentum rapidity y p = tanh −1 (v loc ) smaller than 1.0 for thin (left) and thick (right) shocks. The plots start at the time that a fluid cell attains momentum rapidity 1.0. Clearly, a large fraction of the charge ends up at relatively small rapidities. Due to the hydrodynamic expansion this charge ends up at larger rapidity later on, which explains why the fraction decreases as a function of time.
reaches v z = tanh(1) ≈ 0.76, divided by the total charge. We show this quantity as a function of time in Fig. 10 for both thin and thick shocks.
For thick shocks, the mid-rapidity charge fraction is bigger than for the thin shocks, which implies that the charge distribution is narrower in rapidity, as also illustrated in the right panel of Fig. 4. If the expansion of the system were exactly boost invariant, this quantity would remain constant in the hydrodynamic regime; therefore, the decrease of this fraction at later times is a consequence of non-boost invariant dynamics. Despite this decrease, this computation shows that more than 50% of the charge carried by the initial shocks is deposited within two units of rapidity at early times. 7 This strong baryon stopping is reminiscent of the behaviour of the net proton number in low-energy heavy ion collisions at √ s < 20 GeV. Indeed, heavy ion experiments performed at AGS [18] and SPS [19] found a large stopping of the baryon charge, with most of the net protons concentrated in the mid-rapidity region. As the energy of the beams increases, the rapidity width of the charge distribution increases, reducing the fraction of the total charge at mid-rapidity. Taking the width of our shocks as a proxy for the energy of the incident nuclei [2], this trend is in qualitative agreement with the behaviour of the charge distribution in our simulations. This qualitative agreement indicates that shock wave collisions may provide a good framework to understand the hydrodynamisation of low-and moderateenergy heavy ion collisions, as those studied in the RHIC energy scan. 7 This quantity does not directly correspond to experimental measurements, as this requires a longer and more advanced hydrodynamic evolution and freeze-out. Nevertheless, for comparison, experimentally the fraction of baryon charge between rapidity -1.0 and 1.0 is approximately 37% and 8.4% for heavy ion collisions at √ sNN = 17 and 200 GeV respectively [17].
In contrast, the dynamics of our simple holographic model does not seem to agree qualitatively with the distribution of net baryon number in heavy ion collisions at full RHIC or LHC energies. For those energies, the amount of net baryon number at mid-rapidity represents a small fraction of the total baryon number, and the distribution of baryon charge peaks a few units of rapidity away from the beam rapidity. Quite remarkably, as we have discussed, our thin shock simulations hint at the development of maxima at moderately large rapidity, which may be interpreted as the onset of this non-monotonous behaviour. However, unlike in heavy ion collisions, most of the baryon charge is initially concentrated within a few units of rapidity.
The discrepancy above is consistent with the difficulties in reproducing the LHC multiplicity rapidity profiles at high energies, as noted in [20]. It may be due to the extreme simplicity of our model, which just evolves in the simplest holographic setting with the simplest Maxwell field possible, and for instance does not contain any matter in AdS. We are also working in a model that only describes the very first moments of a high energy collision, without treating a long hydrodynamic phase or freeze-out, and our model restricts to homogeneity in the transverse plane. Naturally, the discrepancy may also point to a deeper difference between our holographic setup and the dynamics of QCD. Indeed, most of the baryonic charge of hadrons are carried by valence quarks, which also carry a large fraction of the full hadron momentum. In contrast, the shock waves posses structure functions concentrated at small Bjorken x [21]. Since processes able to reduce the rapidity of valence quarks by a significant amount involve large momentum transfers, large rapidity shifts at high energies are suppressed as a consequence of asymptotic freedom. The qualitative disagreement in the rapidity distribution of the charge is perhaps not surprising given the absence of this perturbative physics on the gravity side. Nevertheless, since the matter produced at mid-rapidities in a ultra-relativistic collision is soft, the interactions and generation of this matter may still be dominated by strong coupling processes. It would be interesting to develop hybrid approaches able to address these two separated regimes within the duality, perhaps along the lines of [22,23].