Jet Wake from Linearized Hydrodynamics

We explore how to improve the hybrid model description of the particles originating from the wake that a jet produced in a heavy ion collision leaves in the droplet of quark-gluon plasma (QGP) through which it propagates, using linearized hydrodynamics on a background Bjorken flow. Jet energy and momentum loss described by the hybrid model become currents sourcing linearized hydrodynamics. By solving the linearized hydrodynamic equations numerically, we investigate the development of the wake in the dynamically evolving droplet of QGP, study the effect of viscosity, scrutinize energy-momentum conservation, and check the validity of the linear approximation. We find that linearized hydrodynamics works better in the viscous case because diffusive modes damp the energy-momentum perturbation produced by the jet. We calculate the distribution of particles produced from the jet wake by using the Cooper-Frye prescription and find that both the transverse momentum spectrum and the distribution of particles in azimuthal angle are similar in shape in linearized hydrodynamics and in the hybrid model. Their normalizations are different because the momentum-rapidity distribution in the linearized hydrodynamics analysis is more spread out, due to sound modes. Since the Bjorken flow has no transverse expansion, we explore the effect of transverse flow by using local boosts to add it into the Cooper-Frye formula. After including the effects of transverse flow in this way, the transverse momentum spectrum becomes harder: more particles with transverse momenta bigger than $2$ GeV are produced than in the hybrid model. Although we defer implementing this analysis in a jet Monte Carlo, as would be needed to make quantitative comparisons to data, we gain a qualitative sense of how the jet wake may modify jet observables by computing proxies for two example observables: the lost energy recovered in a cone of varying open angle, and the fragmentation function. We find that linearized hydrodynamics with transverse flow effects added improves the description of the jet wake in the hybrid model in just the way that comparison to data indicates is needed. Our study illuminates a path to improving the description of the wake in the hybrid model, highlighting the need to take into account the effects of both transverse flow and the broadening of the energy-momentum perturbation in spacetime rapidity on particle production.

Jets are collimated sprays of particles produced in high energy collisions originating from the production of energetic quarks and gluons which in turn each fragments into a shower of quarks and gluons that subsequently hadronize. These objects play a central role in our understanding of collider physics. In the vacuum, as well as in low density hadronic environments such as pp collisions at the LHC, the physics processes that lead to the formation of these QCD structures are under excellent theoretical and experimental control. This knowledge has transformed those objects into precision tools to understand the physics of the Standard Model and to explore beyond. In dense hadronic environments, like those created in heavy ion collisions at RHIC and the LHC, while our understanding of the properties of these objects are still under intense theoretical investigation , jets are of particular importance since they act as probes of the quark-gluon plasma (QGP) created in those collisions.
The generic phrase used to describe the modification of jets as they traverse a droplet of QGP is "jet quenching". In early studies, the focus was on how the total energy of a jet is degraded; in recent years much attention has been paid to how the shape or structure of the jets are modified in other ways. When jets traverse the plasma, the energetic partons within a jet shower interact with the plasma constituents and, as a result, lose energy to the medium. Many studies have focused on describing how this parton energy loss mechanism takes place. This by itself is a complicated task, since the interaction between energetic partons and the medium that they probe incorporates physics at multiple scales, from the perturbative dynamics of the short distance splitting processes via which the shower of energetic partons develops to the long distance strongly coupled dynamics of the medium.
Since the particular modification of any given jet depends on the parton shower history and varies quite substantially event-by-event, with jets with the same total energy but differing shower structure experiencing very different fractional energy loss and different modifications to their shape, a significant theoretical effort has been devoted to the development of Monte Carlo simulators of jet production events in heavy ion collisions [24,25,[29][30][31].
Complementary to the modification of internal structure of jets, jet quenching also leads to a modification of the medium those objects traverse. This process, which is generically referred to as the backreaction of the medium, is a consequence of the injection of energy and momentum lost by the jet into the plasma. Because momentum is conserved, the wake created in the medium as the jet loses energy and momentum necessarily gives the medium a net momentum in the jet direction, yielding a correlation between the bulk dynamics of the medium and the jet direction. This means that when an experimentalist reconstructs a jet from the hadrons in the final state of the collision, no background subtraction procedure can remove all the particles coming from the hadronization of the medium. What the experimentalist reconstructs as a jet includes hadrons (and energy and momentum) originating from the wake in the medium created by the passing jet as well as hadrons (and energy and momentum) originating from the shower of energetic partons ("the jet itself") after modification by their passage through the medium. Even though the particle production associated with the jet itself (with its modified internal structure) and that associated with the hadronization of the wake in the medium are quite different in character and origin, it is impossible to separate these two contributions to jet observables at colliders. Therefore, a complete description of collider jet data requires an understanding of the dynamics of the backreaction of the medium. Since one of the most striking features of the QGP medium is its almost perfect fluid behavior, a natural expectation is that the dynamics of the "lost" energy, which is to say the dynamics of the wake in the medium, is also well captured by hydrodynamics. These hydrodynamic structures have important consequences for the description of jet observables in heavy ion collisions. Earlier studies applied linearized hydrodynamics with the energy and momentum lost by the jets as external currents, assuming the perturbation caused by the jet energy and momentum deposition into the medium is small [56][57][58][59][60][61]. Many subsequent studies followed this line of thinking [62][63][64][65].
More recent studies relaxed the assumption that the perturbation is linear and carried out more complete analyses using full hydrodynamics [66][67][68][69][70]. A combination of hydrodynamics and weak-coupling partonic transport has also been developed [71]. However, these full hydrodynamic analyses are computationally very expensive.
In this paper we will explore the medium backreaction anew, aiming at a simpler description that can be more efficiently implemented into Monte Carlo generators. While our analysis will not depend on the particular model for jet-medium interactions, we will take inspiration from the hybrid strong/weak coupling model that some of us have developed over the last few years [24,25,29,31,34,36]. This model aims at combining a perturbative description of the parton shower at short distances with a non-perturbative, strongly coupled, description of the interaction between partons in the shower and the strongly coupled hydrodynamic medium via the gauge/gravity duality . In addition, in this hybrid model a very simple prescription for medium backreaction was introduced [29]. This was based on a small perturbation analysis of the medium response on top of an ideal hydrodynamics, to-gether with the assumption that the resulting modification to the final state hadron spectra are small at all p T . Under certain assumptions, which we will review in Section IV, a simple prescription solely controlled by energy-momentum conservation was introduced. Though easy to implement numerically, this prescription for the jet wake particle production in the hybrid model is too soft: comparison to experimental measurements [103][104][105] shows that it yields too many particles with p T < 1 GeV and not enough particles with p T 2 GeV.
(See, for example, Figs. 10, 11, 13 and 14 of Ref. [29].) In this work we will address known oversimplifications of the hybrid model prescription and explore their possible observable consequences. As a first step towards this goal, we apply the framework of linearized hydrodynamics on top of a 1 + 1D boost invariant hydrodynamics (the Bjorken fluid), by assuming the energy and momentum deposited by the high energy partons in the jet shower are small. By solving the linearized hydrodynamic equations numerically, we demonstrate the time evolution of the jet wake and scrutinize the energy and momentum conservation in the evolution. Furthermore, we investigate the effect of viscosity on the bulk dynamics of the jet wake perturbation and on the assumed linear approximation. Then, we convert the jet wake into particle production by using the Cooper-Frye formula. We find that the particle distribution from the linearized hydrodynamics is similar to that from the hybrid model, except that the momentum-rapidity distribution in the former is more spread out. While our hydrodynamic analysis itself does not include transverse expansion, we can nevertheless also explore the effect of transverse flow. We do so by introducing a transverse flow profile, extracted from a 2 + 1D viscous hydrodynamic simulator "VISHNU" [106,107] and calibrated with experimental observables of light particles at low transverse momentum [108], to appropriately boost the jet-induced fluid disturbances and explore the effects of transverse flow on the hadronization of the wake via the Cooper-Frye formula. We find that the transverse flow makes the particle production from the wake harder than in the hybrid model, in line with indications from experimental measurements.
We will also explore how this improvement in the backreaction description could influence jet observables by studying the recovery of energy induced by particles from the wake in conical regions around the direction of the wake (and jet) as well as the momentum distribution of those particles, which affects the fragmentation functions for medium-modified jets.
This paper is organized as follows: In Section II, we introduce our framework of linearized hydrodynamics on a background Bjorken flow. We also explain how to convert the energy loss rate into external currents of energy and momentum in the hydrodynamic equations.
Then in Section III, we explain the numerical recipes to solve the linearized hydrodynamics in detail. In this section we also discuss the jet wake development, tests of energy-momentum conservation, and tests of the validity of the linear approximation. We study distributions of particles produced from the jet wake in Section IV. We also provide comparisons between the results from the linearized hydrodynamics, with and without the effects of transverse flow, and the hybrid model prescription. And, we discuss the dependence of the effects of transverse flow on hadron mass. Next, in Section V we consider the potential effects from the jet wake on jet observables. Finally, we draw conclusions and look ahead in Section VI.

II. LINEARIZED HYDRODYNAMICS ON BJORKEN FLOW
One of the most striking properties of the QGP as produced in high-energy colliders is that it behaves as an almost perfect fluid. This means that the collective properties of bulk particle production, with transverse momentum smaller than a few GeV are correctly captured by a hydrodynamical description of the medium in which all dynamics are reduced to the conservation of the energy momentum tensor, which is approximated by a gradient expansion around equilibrium. However, harder particles may not be assumed to be part of the fluid; as a extreme example, high-energy jets contain fragments of several tens of GeV, which do not have time to arrive close to thermal equilibrium and may be viewed as external probes of the medium. Certainly, the distinction between these two different sets of components is not sharp; a complete description of the jet-medium interaction should be able to account for the dynamics at all energy scales and only below some scale those dynamics would be well approximated by hydrodynamics. However, since currently such description is not available, in this paper we will simply assume that such separation between soft and hard modes exists and that the soft modes are well described via hydrodynamics.
From the point of view expressed above, the energy loss that jets experience in the medium may be viewed as a transfer of energy from hard modes to soft modes. While the overall energy-momentum tensor of the system is conserved, from the point of view of the energy-momentum tensor of soft, bulk particles, T µν this energy-momentum transfer may be viewed as source in the conservation equation, which we may write as where J ν depends on the distribution and dynamics (e.g. energy loss) of hard modes which, in the case of jet-medium interactions, depend on the jet properties.
In the absence of a source, since, as already stated, the dynamics are captured by hydrodynamics, the energy-momentum of soft bulk matter can be written in a gradient expansion.
(See Refs. [109,110] for recent reviews.) The building blocks of this approximation are the (local) energy density field ε and the fluid velocity field u µ , normalized by u 2 ≡ u µ u µ = 1.
The energy-momentum tensor is then constructed using symmetry properties and a powercounting prescription based on the total order of derivatives. To linear order in gradients and neglecting the bulk viscosity for simplicity, we have Here P denotes the equilibrium pressure, which is a function of ε, g µν stands for the metric tensor and η is the shear viscosity. The shear term can be written as where ∆ µν ≡ g µν − u µ u ν .
The main assumption in the hydrodynamic description of the medium backreaction is that the disturbance of the bulk motion of the hydrodynamic fluid (the wake) introduced by the jet can also be described via hydrodynamics. This implies that, also in the presence of the source, the energy-momentum tensor may be approximated in a gradient expansion as in Eq. (2). Certainly, this will be the case sufficiently far away from the jet, when the disturbance has relaxed and is characterised by long-wavelength modes. However, close to the jet the gradients of the flow field will be larger, which may challenge the hydrodynamic assumption. Nevertheless, explicit numerical analysis of the out-of-equilibrium dynamics of both strongly coupled and weakly coupled gauge theories (see the pioneering works [111][112][113][114][115] and see Refs. [92,98,110,116,117] for recent reviews) as well as the recent studies of hydrodynamic attractors [118][119][120][121][122] suggest that hydrodynamics may be applicable even when gradients are large and the matter is away from local thermal equilibrium, although a non-hydrodynamic explanation of the attractor phenomenon in terms of a pre-hydrodynamic epoch during which the dynamics is dominated by one or a few low-lying modes in an effective Hamiltonian would delay the onset of hydrodynamization [123]. Furthermore, explicit analysis of the medium response to an energetic colored object in strongly-coupled gauge theories [79-81, 85-87, 91, 92] demonstrates that, at least in the strongly coupled fluid of N = 4 supersymmetric Yang-Mills theory, the energy-momentum tensor disturbance induced by these energetic particles behave hydrodynamically even at distances as short as 1/πT away from the jet. These theoretical advances motivate us to explore the observable consequences of a hydrodynamic treatment of the backreaction of the medium to high-energy jets in heavy ion collisions. Similar explorations have been performed in Refs. [56][57][58][59][62][63][64][65][66][67][68][69][70].
In addition to assuming the validity of hydrodynamics, in this work we will assume that the disturbance introduced by the external current J µ is small and leads to a perturbation of the stress-energy tensor that can be well-approximated as linear. This assumption allows us to decompose the total stress-energy tensor into two parts: T µν (0) + δT µν , where T µν (0) denotes the stress-energy tensor of the background fluid in the absence of the external current and δT µν is the stress-energy tensor of the perturbation induced by the external current. These two components evolve according to The first of these two equations controls the collective flow of matter in the absence of the jet. The second equation determines the backreaction of the soft modes to the passage of the jet. We want to emphasize that these two equations are not independent from each other. The perturbation piece δT µν generally depends on the details on the background fluid expansion, encoded by T µν (0) . However, the backreaction on the background fluid T µν from the perturbation δT µν can be neglected since the perturbation is assumed to be small.
Because of the interconnection between these two components, to proceed further we need to specify the flow of the background fluid. In a realistic simulation, this bulk flow depends on the particular properties of each event, such as the geometry introduced by the impact parameter of the collision and the fluctuations coming from the nucleonic and subnucleonic structure of the colliding nuclei. While detailed hydrodynamic simulations of the bulk dynamics for such configurations are publicly available, in this exploratory work we will simplify the description of this flow and postpone the analysis of the interplay between the backreaction and the variation of the flow properties in different events to future work.
For this reason, in this work we consider the simplest model for the background expansion: the Bjorken (or boost invariant) flow, which assumes homogeneity in the transverse plane and, therefore, has no transverse expansion. We leave analysis of linearized hydrodynamic perturbations on top of a radially expanding background flow to future studies.
It is easy to define the Bjorken flow in the Milne (τ, x, y, η s ) coordinate system, where η s denotes the spacetime-rapidity and the subscript "s" is used to distinguish it from the shear viscosity. In this coordinate system, the metric tensor is given by g µν = diag(1, −1, −1, −τ 2 ) and the covariant derivative on a (contravariant) vector field can be written as ∇ ν u µ = ∂ ν u µ + Γ µ ρν u ρ . Here ∂ ν denotes the standard derivative and the nonvanishing Christoffel symbols are given by The velocity field in the Bjorken flow is given in these coordinates by Then one can work out the shear term in the stress-energy tensor: 2∇ µ u ν 0 = diag(0, 2 3τ , 2 3τ , −4 3τ 3 ). The hydrodynamic equation ∇ µ T µν (0) = 0 for the background Bjorken flow is given by To solve the hydrodynamic equation (8), we need the equation of state that relates P 0 to ε 0 . We will use the textbook result for a non-interacting gas of gluons and N f = 2 flavors of massless quarks, which is not far from the results of lattice QCD calculations of the QCD equation of state for QGP with N f = 2 + 1 (2 light and 1 strange quarks) at high temperature [124]. For the background energy density and pressure, we have where s 0 is the entropy density and g = 40 is the number of massless degrees of freedom in non-interacting N f = 2 QGP. This equation of state leads to a speed of sound which also controls the dynamics of the fluid perturbation.
To solve the system, we also need to specify the value of the viscosity. In this paper, we will consider and compare two characteristic values of this transport coefficient: the ideal case, η = 0; and the value inferred from strongly coupled computations η = s 4π , where s denotes the entropy density.
To explicitly write out the hydrodynamic equation ∇ µ δT µν = J ν for the perturbation, we must work out the explicit expressions for the perturbed stress-energy tensor δT µν and the external current J ν , which is what we are going to do in the next two subsections.

A. Linear Perturbation of Bjorken Flow
As discussed earlier, we assume the perturbation caused by the jet energy and momentum deposition leads to a small perturbation. So, the perturbed velocity field can be written as Since the perturbation δu µ is small, the fluid velocity is still properly normalized if we neglect higher order terms beyond the linear order: (u µ 0 + δu µ ) 2 = 1 + O((δu) 2 ).
The total stress-energy tensor T µν with the perturbation can be calculated by plugging u µ 0 + δu µ into Eq. (2) and decomposing ε = ε 0 + δε and P = P 0 + δP . Since ε + P = T s and η = s 4π , we also decompose the shear viscosity where δε = T δs and we have defined The physical significance of γ η is that it controls the relaxation time of hydrodynamic excitations of the fluid. Expanding the expression (2) for T µν to linear order in the perturbation and subtracting the unperturbed background T µν (0) , we find where We use the short-hand notation ∂ x = ∂ ∂x and similarly for y, η s and τ . The linearized hydrodynamic equation ∇ µ δT µν = ∂ µ δT µν + Γ µ ρµ δT ρν + Γ ν ρµ δT µρ = J ν can then be written as where the bold symbols denote Euclidean 2-vectors and the dot product between two Euclidean vectors is given by a · b = a x b x + a y b y .
We solve the equations of linearized hydrodynamics in momentum space. For simplicity, we define The Fourier transformed quantities in momentum space, which we denote by δε,g ⊥ ,g η , and J µ are then defined via where we are slightly abusing the notation: we use g η instead of g ηs in the definitions and the momentum conjugate to the spacetime-rapidity η s is labeled k η rather than k ηs . We believe the physical meanings of g η and k η are clear; they have nothing to do with the shear viscosity η.
Plugging the Fourier transforms into the linearized hydrodynamic equations and using the definition of the speed of sound, we obtain the equations for the excitation in Fourier space. In these equations, we will perform one additional approximation: we neglect all terms that are suppressed by a factor of order γ η /τ relative to terms that are of order unity, with γ η as defined in (13). This approximation is motivated by the origin of the hydrodynamic approximation as a gradient expansion and by observing that due to the boost invariant symmetry assumed in Bjorken flow, the gradients of the velocity field scale as 1/τ . In essence this approximation amounts to considering the flow field in regions where γη τ 1, where the strict gradient expansion is valid. With this approximation, the equations for the perturbations are For later convenience, we decomposeg ⊥ into the longitudinalg L and transverseg T modes: Then from (28) we obtain the hydrodynamic equations forg L andg T : These equations illustrate the convenience of the decomposition of the momentum flux we have just introduced, since we see that the hydrodynamic equation for the transverse mode decouples from the other equations. This is the manifestation of the diffusive mode of a homogeneous plasma in this expanding system. In contrast, the longitudinal mode couples with the energy perturbation δε and the η s -momentum perturbationg η . This coupled channel describes the sound excitations of the system. These two channels are sourced by different components of the current J µ . To be able to solve the dynamics, we need to specify the functional form of this source.

B. External Current from Hybrid Model
In this subsection we describe how we treat the energy deposition by the jet into the plasma and how we use this knowledge to model the source term that feeds the hydrodynamic response of the medium. While detailed calculations of the energy loss rate of hard partons propagating plasma are available in the literature, the hydrodynamization processes via which this energy becomes a modification to (e.g. a wake in) the hydrodynamic fluid and the manifestation of these deposition and hydrodynamization processes via this type of source terms is not under good theoretical control. For this reason, the specification of this source will require additional assumptions. Whenever possible, we will draw inspiration from explicit computations of the jet-medium interaction in both the weak and strong coupling limits.
Although some of our consideration will be general, to make the analysis concrete we will focus on the strong-coupling limit of jet-medium interaction and use knowledge from the N = 4 SYM description via a gravity dual. Following the hybrid strong/weak coupling model perspective, we will employ the following energy loss rate for a single hard particle traversing the fluid [97,99] where E in is the "incoming" energy of the parton, before it loses any energy to the hydrodynamic medium, and T is the local temperature of the QGP medium. In strongly coupled N = 4 SYM, this energy loss rate is fixed and the parameter κ sc = 1.05λ 1/6 , with λ the 't Hooft coupling. In the hybrid strong/weak coupling model it is assumed [24] that the functional form of this rate remains the same in strongly coupled QCD as in strongly coupled N = 4 SYM gauge theory and that all differences between these two theories can be condensed into the value of this parameter, which becomes a fitting parameter of the model. Taking the phenomenological analysis of Ref. [34] as a guide, we shall choose the value κ sc = 0.4 throughout this paper.
In a future Monte Carlo study, we would apply the rate of energy loss (34) and the entire analysis of the present paper to each energetic parton in each jet shower in an ensemble of events in which in each event jets are produced at a different location and propagate in different directions (typically two or more jets, which may have interesting consequences [65,125]) within an expanding, cooling, droplet of QGP. In the present study, we shall simplify things by considering only a single energetic parton traversing the (boost invariant; Bjorken) hydrodynamic fluid.
To fix the total energy-momentum injection from the hard parton into the plasma we also need to determine the rate of momentum loss. We can easily relate this rate to the rate of energy loss by assuming that the parton moves at the speed of light with a fixed spacetime rapidity η s parton and transverse directionn ⊥ parton , These rates then constrain the functional form of the current J µ that describes the injection of energy and momentum into the medium but do not by themselves fix it, as we shall see.
Working in Milne coordinates and assuming that the disturbances induced by the jet fall sufficiently fast at spacetime rapidities away from the jet, we can begin by relating J µ to the rates of energy and momentum loss by integrating over the space-like fixed-τ hypersurface.
By means of Stokes' theorem The integral constraints (39) and (40) will serve to fix the magnitude of J µ but they do not fix the full functional form of the current. Following physical reasoning, the simplest assumption is that most of the energy deposition occurs in some small region around the local position of the hard parton while it is traversing the plasma. Assuming a Gaussian parametrization for the deposition process of energy and momentum in the rest frame of the fluid cell, we write the current (in Milne coordinates) as where (x parton (τ ), y parton (τ ), η s parton ) is the transverse position at time τ , and the spacetime rapidity, of the energetic parton propagating through -and losing energy and momentum to -the medium, and where σ and σ ηs are the Gaussian widths in the transverse and spacetime-rapidity directions, respectively. With this parametrization and for fixed values of the width parameters, the τ -dependent functions C µ (τ ) are fixed by the rates of energy and momentum loss through the integral constraints (39) and (40). Using the relations between the Minkowski t and z components of J µ and the Milne τ and η s components, as well as Eqs. (39,40), we obtain After fixing these coefficients we have fully specified the source term and we are now ready to solve the linearized hydrodynamics equations, (27,29,32,33,34,41), to study the jet wake. In the next Section, we will explain the numerical procedure that we use to solve these equations and to generate the medium response to a high-energy parton.

A. Numerical Setup
In the previous Section, we have completely specified the equilibrium and transport properties of the system we simulate as well as the rates of energy and momentum loss that govern the energy and momentum injected into the medium and the equations that govern the medium response. Prior to solving these equations and describing the backreaction of the medium numerically, we also need to describe the geometry and parameters that specify the events we simulate.
As we have said, in this paper we shall choose a particularly simple hydrodynamic solution for the bulk fluid, treating it as a longitudinally expanding boost invariant fluid with no transverse expansion. We take the initial temperature of the fluid to be T 0 = 400 MeV at an initial time that we choose to be τ 0 = 0.6 fm/c. The bulk matter then evolves hydrodynamically, expanding longitudinally and cooling, in a way that is determined by the equation of state and shear viscosities specified above. We shall assume that the background flow freezes out, which is to say hadronizes, when the local temperature drops to a freezeout temperature that we choose to be T f = 154 MeV. The freezeout time is τ f ≈ 10.5 fm/c for the ideal case in which we choose η = 0 and is τ f ≈ 11.7 fm/c when we choose a fluid with nonzero shear viscosity, specifically η 0 = s 0 /(4π). When the fluid viscosity is nonzero, heat can be generated as the fluid expands; this is why it takes a little longer for the fluid to cool to T f . While the first order hydrodynamic equation (8) is known to lead to unphysical behavior due to the acausal nature of the relativistic Navier-Stokes equations, for the values of η 0 /s 0 , τ 0 and T 0 that we use in this work -which satisfy τ 0 T 0 η 0 /s 0 -these effects are negligible and the solution of the equation (8) provides a good description of the evolution of the system obtained with causal completions of hydrodynamics [126].
On top of this fluid, we now wish to compute the linearized hydrodynamic response to an energetic probe traversing the medium and losing energy at the rate (34), with the strong coupling parameter taken to be κ sc = 0.4. We shall choose the Gaussian widths in Eq. (41) as σ = 1 πT and σ ηs = 1 π . We will study a high energy parton with an initial energy E in = 100 GeV that starts at a point (x 0 , y 0 ) in the transverse plane and travels along the +x-direction in the transverse plane and η s parton = 0. We assume the jet direction does not change during the in-medium evolution, so x parton (τ ) = x 0 + 1 cosh(ηs parton) (τ − 0.6) fm, y parton (τ ) = y 0 . The shape of the solution is independent of the choice of x 0 and y 0 , since our fluid is translationinvariant in the transverse directions. In a future Monte Carlo study we would choose the probability distribution for x 0 and y 0 for an ensemble of jets using a Glauber model and would choose their direction of propagation at random. In the present calculation, with a medium that is homogeneous in the transverse directions, any choice of x 0 and y 0 is as good as any other, as long as we make sure that the wake generated is contained within the box that we use for our numerical computation.
Later, in Section IV, we shall model the effects of transverse, radial, flow by using the velocity field extracted from a longitudinally boost invariant 2 + 1D viscous hydrodynamic simulation of a central heavy ion collision to boost the hadrons that freezeout from each point on a freezeout surface. Doing so will mean that each point in the medium is no longer equivalent to others, and in particular will require choosing a point in the transverse plane that is the origin for the radial flow, which need not be the same as the point where the high energy parton starts. This also means that, in Section IV, it will be necessary for us to make sure that (most of) our wake falls within the transverse extent of the freezeout surface extracted from the 2 + 1D hydrodynamic simulation.
We introduce an additional parameter τ e , after which the jet energy loss is switched off.
This parameter is necessary since our system is infinite in extent and homogeneous in the transverse plane and, unlike in a real collision, jets never leave the plasma. To mimic the finite size of a real collision, we impose that for τ > τ e , the right hand sides of (27,28,29) vanish. We can change the value of τ e to change how much energy the jet loses during the in-medium evolution.
Once all parameters are specified we are ready to simulate the medium response numerically. We do the calculation in momentum space, solving the Fourier transformed versions of the linearized hydrodynamics equations, namely Eqs. (27,28,29). We discretize momentum In practice, we only need to solve the linearized hydrodynamics equations for points in mo- 20,20] since it turns out that contributions from |k η | ∈ [20,60] are already exponentially suppressed.
Our choice of the momentum range for k x gives a resolution of ∆x = π (k x )max ≈ 0.06 fm in the transverse plane with (k x ) max = 10 GeV, and similarly for the y-direction. For the η s -direction, the resolution power is ∆η s = π 60 ≈ 0.05. The range of the x, y, η s space included in our calculation is a box of size πNx (k x )max × πNy (k y )max × πNη (k η )max centered at the origin. That is, we establish as our convention that x 0 = y 0 = 0 corresponds to the center of our calculational box. We must then choose N x , N y and N η large enough such that our calculational box encompasses the whole jet wake in the transverse plane and spacetime rapidity. The conditions can be approximately written as where the left-hand-side of the first line gives the farthest distance the jet wake can reach along +x direction while left-hand-side of the second line gives that along −x direction. In the viscous case, the size of the jet wake can be somewhat bigger than the above estimate because perturbations to the hydrodynamic fluid can diffuse in addition to spreading at the speed of sound. Since the farthest distances along the +x and −x directions are different, we can optimize the value of x 0 to minimize the required value for N x , which further reduces the computing resources needed to solve the linearized hydrodynamics. After fixing N x , we choose N y = N x by symmetry.
We shall present results from two calculations in this paper: temperature of the background medium has dropped to T f = 154 MeV, and the second time plotted is at the midpoint in time between the first and third. The wake moves along the +x-direction, with energy accumulation on a wavefront that is apparent in the Figures and that follows the energetic parton, and energy depletion just behind this wavefront. In the case of an ideal fluid, Case 1, only the propagating sound mode exists while in the viscous fluid, Case 2, there is also a diffusive mode. This is why the energy perturbation in the viscous fluid, Case 2, spreads out more in the transverse plane, with its peak value dropping more rapidly with time. Note that at each time plotted the left and right panels in Fig. 1 have the same vertical axes, but the vertical axes for plots at different times differ. We see that the magnitude of the perturbation is initially similar in the two cases, but indeed it drops much more rapidly with time in the viscous fluid.
Next we plot the momentum perturbation g x (x, y, η s = 0) in Fig. 2. Since the jet moves along the +x-direction and only deposits momentum along the +x-direction, we will not show plots for g y here. In the ideal case, a wing-shaped peak structure in g x is left behind by the energetic parton propagating along the x-axis, and this moving wake of fluid remains quite prominent at late times. This structure arises from the transverse mode g T . The transverse mode does not diffuse in the ideal case (see Eq. (33)). In addition to sustaining the wing-shaped peak along the x-axis, a Mach cone structure is also seen in a wavefront spreading out from the place where the energy density has its positive peak in Fig. 1.
This structure originates from the coupling between the longitudinal mode and the energy perturbation (see Eq. (32)). In the viscous case, the transverse mode g T diffuses and spreads so we see much less accumulation of momentum perturbation along the x-axis, a much less prominent wing-shaped structure along the +x-direction, less wake. For the same reason, the Mach cone structure coming from the longitudinal mode is also spread out. Note again that at each time plotted the left and right panels of Fig. 2 have the same vertical axes (whereas we have chosen different vertical axes at different times). We see that, as in Fig. 1, the magnitude of the perturbation drops much more rapidly with time in the viscous fluid.

C. Energy-Momentum Conservation
As a test of our numerical solver, we check that the energy and momentum injected into soft modes by the passing high energy parton are indeed stored in the hydrodynamic fields at the (much later) time of freezeout. As we discussed in Section II B and Appendix A, the integration of δT τ τ = δε and δT ηsτ = g η at some τ > τ e do not reproduce the energy and the z-momentum deposited by the high energy parton. Instead, the conserved energy and momentum are related to the integrals of the stress tensor components as τ dx dy dη s δT tτ (τ, x, y, η s ) = ∆E tot (51) τ dx dy dη s δT ⊥τ (τ, x, y, η s ) = ∆P ⊥ tot (52) τ dx dy dη s δT zτ (τ, x, y, η s ) = ∆P z tot , where δT tτ = (cosh η s )δT τ τ + τ (sinh η s )δT ηsτ and δT zτ = (sinh η s )δT τ τ + τ (cosh η s )δT ηsτ .
In To complete our check of energy and momentum conservation, we calculate the total energy and momentum lost by the high energy parton from the strong coupling formula (34) together with (36): in which the velocity v i is given in Eqs. (37,38). The quantities ∆E tot and P i tot represent the total energy and momentum lost by, and hence deposited into the medium by, the passing high energy parton.
For each of the two Cases considered here, we calculate the total deposition of energy and momentum into the fluid by the energetic parton and compute the energy and momentum in the hydrodynamic wake at freezeout from the linearized hydrodynamics. Since in both our calculations we have a high energy parton moving along the +x-direction with constant    η s = 0, we have v x = 1, v y = 0 and v z = 0. Therefore, all the momentum deposited is along the +x-direction. In Table I we show our results for the energy and momentum in the x-direction stored in the hydrodynamic fields at freezeout and compare them with the energy and momentum deposited into the plasma by the high energy parton. We observe that indeed energy and momentum conservation both work well. We also observe that momentum conservation is more accurate than energy conservation. We believe this is due to the cosh η s and sinh η s terms in T tτ , which can lead to numerical instability at large η s .
Our test of energy conservation confirms that we have this well under control. We have also tested that momenta along the two directions perpendicular to the jet are very small in magnitude (< 10 −15 ), consistent with zero, and very well conserved.

D. Validity of Linear Approximation
To close this Section, we check the validity of our use of linearized hydrodynamics, which is based on the assumption that δε/ε 0 (and hence δP/P 0 ), δu x , δu y and δu ηs are all small, meaning that a truncation to linear order in these quantities is a good approximation.
To this end, we show the contour plots of δε/ε 0 and δu x at different times in Figs These observations also highlight the importance of incorporating viscous effects into the dynamics of the wake, since turning on a viscosity as small as that in Case 2 introduces quite substantial changes to the magnitude of the wake by the time of freezeout. We have reported results for Case 1, the ideal fluid, because it is a standard calculational benchmark.
In the next Section, however, when we consider particle production from the perturbed fluid we shall focus entirely on Case 2, where the fluid has a nonzero viscosity and the linearized approximation that we employ throughout is under good control long before freezeout.
For later reference, we plot δu y (η s = 0) and δu ηs at freezeout in Fig. 6. Since δu ηs = 0 at zero rapidity because of symmetry, we show its value at η s ≈ 0.52. As those plots clearly show, δu y and δu ηs are much smaller than δε/ε 0 and δu x , consistent with the discussion above.
Finally, in Fig. 7 we plot the vector field δu ⊥ = (δu x , δu y ) at freezeout in the transverse plane at η s = 0. We see a vortical flow pattern in the wake left behind by the energetic parton in both the ideal and viscous cases. As discussed in Ref. [127], the jet creates a vortex freezeout hypersurface that we shall denote by σ µ yields the spectrum for hadrons with a given mass m produced at freezeout: where f (x) is the Bose-Einstein distribution (e x − 1) −1 for bosons with mass m and the Fermi-Dirac distribution (e x + 1) −1 for fermions with mass m.
The freezeout hypersurface is typically taken to be the isothermal hypersurface at which the fluid temperature drops to a specified T f . Within our assumption of homogeneous Bjorken flow, for the unperturbed background fluid the isothermal freezeout surface is also isochronous, since all fluid cells cool down to the freezeout temperature T = T f at the same freezeout time τ f . For the bulk particles, taking the velocity field as u µ 0 = (u τ 0 , u x 0 , u y 0 , u ηs 0 ) = (1, 0, 0, 0), we find where m T ≡ p 2 T + m 2 is the transverse mass (p T is the transverse momentum and m denotes the mass of a specific particle produced from the wake) and y is the momentum rapidity.
When we consider the perturbed fluid, including the wake that originates from the deposition of energy and momentum by the passing high energy parton and that has evolved hydrodynamically until the freezeout time as we have analyzed in previous Sections, the isochronous and isothermal surfaces no longer coincide because the energy perturbation δε changes the value of the temperature at the freezeout time of the background fluid. Nevertheless, in the case of a viscous fluid (even with a viscosity as small as η 0 /s 0 = 1/(4π) as in Case 2), at the freezeout time τ f the perturbations from the wake are small, which means that the consequent change in the temperature of the fluid in the wake relative to the background temperature is also small and so are the differences between the isochronous and isothermal surfaces. For this reason, throughout our analysis we shall compute particle production at the fixed τ f , isochronous, hypersurface. The perturbations to the fluid at the τ = τ f hypersurface due to the presence of the wake perturb the final-state spectra of hadrons, as we shall compute.
At the time τ = τ f , the flow velocity is u µ = (1, δu x , δu y , δu ηs ), where the velocity perturbations are related to g ⊥ and g η via Eqs. (21,22), and g ⊥ and g η can be extracted from our linearized hydrodynamic calculation at τ = τ f . Also, the local temperature at freezeout is modified by δε, see Eq. (9). From Eqs. (9,10,21,22) we then find where g = 40 is the number of degrees of freedom in the plasma phase. By taking the difference between the particle distribution from the fluid with the perturbation and that without, we obtain the distribution of particles originating from the jet wake in the linearized hydrodynamic approach: where φ is the azimuthal angle of the transverse momentum of the produced particle. The determination of the spectrum and the distribution in azimuthal angle and momentum rapidity y from the expression (61) requires knowledge of the hydrodynamic perturbations, which in the linearized approximation to hydrodynamics we obtain via the numerical results from the previous Section.
In the hybrid strong/weak coupling model of Ref. [29], a series of further assumptions and simplifications were employed, which reduced (61) to a closed-form expression depending only on the amount of energy and momentum deposited in the fluid and not on knowledge of the detailed form of the resulting fluid perturbation. The analysis of Ref. [29] was, first of all, based on ideal hydrodynamics where the comoving entropy τ s (s denotes the entropy density) is conserved. Second, in the Cooper-Frye analysis of freezeout (56), the distinction between fermions and bosons was ignored in Ref. [29], with the distribution f taken as the Boltzmann distribution. Third, in Ref. [29] the Boltzmann exponential was expanded to linear order in the velocity perturbation, which is to say that it was assumed that small perturbations to the hydrodynamic fluid yield small perturbations to the hadron spectra at all p T , an assumption that is only valid at low p T . We shall not need to make any of these simplifying assumptions in our treatment of the wake. The fourth assumption made in Ref. [29], which we shall also assume, is that the direction of the high energy parton, in azimuth and in spacetime rapidity, does not change as the high energy parton propagates through the fluid. In particular, this means that the high energy parton does not deposit any momentum in the η s direction into the fluid, see (46). Fifth, it was also assumed in Ref. [29] that the energy and transverse momentum lost by a jet with momentum rapidity y J turns into perturbations of the fluid at spacetime rapidity η s = y J throughout, i.e., δε(τ, η s ) ∝ δ(η s − y J ) and δP ⊥ (τ, η s ) ∝ δ(η s − y J ) for all τ , including at τ = τ f , the freezeout time. In our linearized hydrodynamics calculation, we shall not need this assumption either.
Making these five assumptions and simplifications, as in Ref. [29], yields the following hybrid model expression for the excess of particles generated by the wake in the fluid at freezeout: d∆N hybrid p T dp T dφ dy = 1 32π where y J and φ J are the momentum rapidity and transverse azimuthal angle of the "jet", e.g. the high energy parton.
In this work we will test how particle production described by the oversimplified hybrid model expression (62) compares to particle production described by the full linearized  56) and (61), this effect on the particle spectrum can be very significant and is certainly much more important than the first effect above. In this Section, we will estimate the modification of the spectrum and the azimuthal and momentum rapidity distributions of the hadrons originating at freezeout due to this second effect.
According to the previous discussion, in our estimate of the effect of transverse flow we will assume that the flow profiles that we have computed in the previous Section remain unchanged, as does the freezeout hypersurface. However, we will imagine that these profiles, which is to say the wake that we have computed, freezes out from a hypersurface at which the radial, transverse flow is significant, meaning that when the particle production needs to be calculated in the local rest frame the boost to this frame now includes a transverse boost.
When the fluid breaks into particles, the momentum of a particle produced at each cell in the local rest frame, p µ cell , whose distribution is determined via the integrand of equation (56), differs from the momentum observed in the laboratory frame, p ν lab by a local transverse boost where Λ µ ν (v cell ) is a Lorentz transformation with v cell = v cell (τ, x, y, η s ) which is the transverse velocity of the fluid cell seen in the laboratory frame. The longitudinal boost has already been accounted for when we write u 0 = (1, 0, 0, 0) in the (τ, x, y, η s ) coordinate system, so we only need to include the transverse boost in Eq. (63). This modification is incorporated into the distribution of particles produced from a given fluid cell at freezeout by replacing the u µ 0 p µ and u µ p µ terms in the distribution (61), respectively, by their boosted counterparts: We then need to obtain v cell for each fluid cell on the freezeout surface from a longitudinally boost invariant 2 + 1D hydrodynamic calculation.
For our estimates, we extract the spacetime dependence of the transverse flow from a hydrodynamical simulation that incorporates radial transverse flow obtained using the "VISHNU" package [106,107], which solves longitudinally boost invariant 2 + 1D viscous hydrodynamics for azimuthally symmetric, central, heavy ion collisions. We shall use results from the calculations calibrated with experimental measurements on low transverse momentum observables found in Ref. [108]. The extracted transverse velocity at mid-rapidity v cell = (v x , v y , 0) at the freezeout hypersurface determined by T f = 154 MeV is depicted in Fig. 8 for central collisions (centrality 0-5%). When the transverse flow is taken into account, the freezeout hypersurface is finite in extent in the transverse plane. Outside its boundary, we simply set the transverse flow velocity to be zero. To make sure that this oversimplification has little influence on our results, we must ensure that the wake in the hydrodynamic fluid is located, at the time of freezeout, almost entirely within the colored circles in Fig. 8. We choose the centers of these circles, e.g. the point of origin of the radial flow at the center of the collision, to lie at the center of our calculational box. And, the choices that we made in Section III for the point of origin of the high energy parton and the direction of its motion, relative to our calculational box, were made with malice aforethought: they ensure that the wake that we have calculated lies within the colored circles in Fig. 8. (See the bottom-right panels in Figs. 4 and 5 and the right panels of Fig. 6.)

C. Results, Including the Effects of Radial Flow
We are now ready to look at the results of our calculations of the distribution of particles coming from the freezeout of the perturbed hydrodynamic fluid. We only consider what we referred to as Case 2 in Section III, which is to say the case where the perturbation is the wake that a high energy parton leaves in a fluid with a small but nonzero viscosity.
For the present, we assume that all the particles produced at freezeout are pions with mass m π = 138 MeV. We begin by looking at the spectrum of pions produced from the wake, shown in Fig. 9. Then, in Fig. 10  ∆E tot = dp T dφ dy m T cosh y d∆N dp T dφ dy .
First, we compare the p T spectra shown in Fig. 9. The linearized hydrodynamic calculation with no effects of transverse flow yields a similar p T spectrum as that obtained in the simplified hybrid model calculation, although the normalization is different because the momentum rapidity distribution in our linearized hydrodynamic calculation is significantly broader (see below). The similar p T shape is expected since the particle spectrum in both cases must be close to a thermal spectrum at T f = 154 MeV when the perturbation is small, as is the case when the fluid has a nonzero viscosity which allows for the perturbation to spread diffusively, as we have seen in Section III D, We also see in Fig. 9 that the inclusion of radial transverse flow leads to a blue-shift in the effective temperature, which makes the p T spectrum harder, as compared with the other two cases. In addition, we notice that the modifications introduced by radial flow make the spectrum of the wake negative at low p T . What we are computing -see (61) -is the difference between two spectra: one for pions produced at freezeout from the fluid perturbed by the wake, and the other for pions produced at freezeout from the unperturbed fluid. In kinematic regions where this quantity is negative, there are fewer pions after we include the perturbation induced by the wake left behind by the high energy parton. What we see from Fig. 9 is that this does not arise without radial flow, but the boost created by radial flow hardens the spectrum of the wake so much that in addition to pushing it up significantly at p T 1 GeV it makes it negative at the lowest p T .
Next we analyze the distribution in the azimuthal angle, shown in the upper panels of The shape of the inclusive azimuthal distribution remains also similar when transverse flow effects are included. However, since the p T -distribution is severely modified, it is interesting to compare the different distributions in azimuthal angle for semi-hard, p T > 1 GeV, particles as this is the kinematic regime in which radial flow acts to push the spectrum up (see we include all pions whereas in the right panels we include only those pions that have see this, consider Eqs. (64) and (65), substituted into (61).) In a future realistic simulation, however, jets would be produced at all starting points and in all directions, and although the parton direction and the flow direction will tend to be aligned by the time of freezeout there will certainly be jets in the ensemble for which this is not the case. Therefore, when all such configurations are averaged over in a future more realistic Monte Carlo study, the magnitude of the observed correlation (e.g. the degree to which the azimuthal distribution becomes more peaked due to the effects of radial flow) is likely to be somewhat reduced.
Next, we discuss our results for the distribution of particles from the wake in momentum rapidity, depicted in the second row of Fig. 10. For this distribution, the linearized hydrodynamic calculation with no effects of radial flow and the hybrid model approximation differ significantly. The linearized hydrodynamic calculation, without transverse flow, leads to a much broader distribution in y than the hybrid model. This originates from the broad distribution of the energy and momentum of the wake in spacetime rapidity η s in our linearized hydrodynamic calculation. Once the energy and momentum are deposited, they propagate along the ±η s directions due to the excitation of sound modes (see Fig. 3). This is in contrast to the hybrid model assumption, that considers the energy and momentum perturbation is located at a fixed η s throughout. This difference in rapidity distribution is also behind the difference in yield between these two sets of calculations. Since we normalize the spectrum by conserving the total energy (66), which involves the integration of m T cosh y, the height of the p T spectrum from the linearized hydrodynamic calculation without any effects of flow is smaller than that from the hybrid model, see Fig. 9, because of the broader momentum rapidity distribution in the former, see Fig. 10. However, while the distribution obtained from the linearized analysis is rather wide in the absence of flow, incorporating the effects of radial flow reduces its rapidity extent, as particle production from the wake becomes more strongly correlated with the direction of motion of the high energy parton.
Finally, we note by comparing panels (b) and (d) of Fig. 10 that introducing radial flow substantially increases the number of particles produced with p T > 1 GeV, confirming what we already saw in Fig. 9.
To conclude this Section, we explore a natural consequence of transverse flow. In Fig. 11, hadron increases. This is a standard tell for hydrodynamic behavior. We see in Fig. 11 that the same phenomenon occurs for the particles coming from the wake. For hadrons with a higher mass, the location of the positive peak in the p T -spectrum is pushed to a higher value, which is expected since the momentum gained by the particle from the transverse flow grows with the particle mass. This is not a surprise since we are treating the wake via linearized hydrodynamics, but it is good to see. This may lead to interesting observables in the fragmentation function. For example, the ratio of the pion fragmentation functions in pp and AA collisions may have a different shape from that of the proton fragmentation functions, since more protons at 2 GeV are produced than pions at 2 GeV from the jet wake.

V. TOWARDS JET OBSERVABLES
As has been already demonstrated by several theoretical analyses [26,29,35,36,71,125,[129][130][131][132][133][134][135][136][137], some jet observables can be significantly modified by medium response effects. By momentum conservation, the wake that a jet leaves behind in the droplet of QGP through which it has passed must carry momentum in the jet direction, since it carries the momentum lost by the jet. That means that after hadronization the particles originating from the wake also must have net momentum in the jet direction. And in turn this means that when experimentalists define what they observe as a jet -via a reconstruction algorithm and a background subtraction procedure -the jets that they reconstruct must include some particles coming from the wake since background subtraction procedures can at best only subtract effects of fluctuations that are uncorrelated with the jet direction. These effects [138,139] are likely comparable to or larger in magnitude than the effects originating from the wake of the jet, highlighting the importance of background subtraction in any experimental analysis. In an experimentally reconstructed jet, the harder particles likely originate directly from the jet shower whereas some of the soft particles originate from the wake in the plasma.
So, jet observables that are sensitive to the softer particles within jets are likely to be significantly influenced by the wake in the plasma.
The contribution of particles in jets that originate from the wake in the plasma to jet observables has been particularly well studied in the hybrid strong/weak coupling model.
In this model, the inclusion of these particles coming from the backreaction of the medium is essential to the description of many basic jet observables, like jet suppression and the medium-modification of fragmentation functions [29]. The calculations of the previous Section, in particular those illustrated in Figs. 9 and 10 (b), show that with our better description of medium response built upon linearized hydrodynamics and incorporating the effects of radial flow, when the wake of a jet hadronizes it yields more semi-hard particles (1 GeV p T 3 GeV) that are more correlated along the jet direction than in the simplified implementation of the wake used in the hybrid model. At a qualitative level, both these features are exactly what is indicated by the discrepancies between the hybrid model calculations of so-called missing-p T observables in Ref. [29] and the experimental measurements of these observables in heavy ion collisions at the LHC by CMS [105]. This provides strong motivation for a future quantitative study of the consequences of this improved description of jet wakes for jet measurements at the LHC. To make quantitative comparison with experiments, however, this improved description of jet wakes will need to be incorporated into a Monte Carlo study of a realistic ensemble of jet showers with a distribution of points of origin, direction, and shower structure, each propagating through a realistic model for the droplet of QGP, for example as in the hybrid model. We leave this to future work.
Nevertheless, to gauge the influence that this improved description could have on actual jet observables, in this work we will study two examples which share features with some of those jet observables. In each case, we will compute only the contribution to the observable that comes from the hadronization of the jet wake sourced by our single high energy parton, as in Section IV. In particular, we will compute the energy carried by particles coming from the backreaction of the medium that lie within a cone of a specified opening angle centered on the direction of the high energy parton; and the energy distribution of those particles.
For simplicity, as in Section IV we will only consider a setup in which the high energy parton that loses energy propagates at zero rapidity; we shall measure azimuthal angles with respect to its transverse direction.

A. Energy Recovered inside a Cone
Our first example is the amount of energy lost by the high energy parton carried by particles originating from the wake in the plasma recovered within a cone of opening angle dφ dy dp T p 2 T + m 2 cosh y d∆N dp T dφ dy .
In any realistic simulation, any jet reconstruction algorithm will collect particles around the jet direction and incorporate them into the definition of the jet. If backreaction occurs, particles from the medium response will also be added to the jet. As a consequence, even if there is a transfer of energy from hard jet fragments to soft medium particles, those medium response particles will put energy back into the reconstructed jet. Therefore, Eq. (67) gives us a way to gauge how much of the "lost" energy is "recovered", namely included as a part of the reconstructed jet coming from the backreaction of the medium. Our results for this quantity as a function of R are shown in Fig. 12. The plots extend out to large enough values of R 2 that the total amount of energy lost by the high energy parton is fully recovered.
We note from the solid magenta curve in Fig. 12 that while the hydrodynamic calculations recover the full deposited energy, a very small fraction of the energy is missed in the hybrid model calculation, even at very large angles. This is just a consequence of the approximations used to derive the simplified expression, Eq. (62), which assumes the production of massless final particles, but is used to study the productions of pions.
To gauge the importance of semi-hard medium response particles in ∆E(R), in Fig. 12 we also show results when we restrict the transverse momentum of the recovered particles to p T > 1 GeV. Doing so highlights an important comparison, namely the very different momentum composition of the medium response particles within a given angular region. For both the calculations with no transverse flow, energy is mostly recovered in soft particles, and the contribution from semi-hard particles, p T > 1, to ∆E(R) in Fig. 12 is very small.
However, when the effects of flow are added this distribution completely changes and a large fraction of the energy within the cone is recovered in the form of semi-hard particles, which   is consistent with what we saw in the spectra of Fig. 9. We will further investigate this in the next Subsection.
First, though, we should understand one of the most striking features of the energy recovery plots in Fig. 12, namely the non-monotonic behavior observed in all three calculations.
The fact that ∆E(R) is above its large-R value for cones with opening angle R ∼ 2 means that after hadronization the wake deposits more than its total energy, more energy than was lost by the high energy parton, within a cone with this radius. How can this be? This is a consequence of the fact that the spectrum in Fig. 9 and the azimuthal and spacetime rapidity distributions in Fig. 10 are all negative in some kinematic regions. Recall that this means that in such regions the presence of the wake depletes the energy/momentum distribution relative to what it would be in the absence of the wake, which we can understand as a result of a physical effect: the high energy parton drags the fluid along its direction of propagation, creating a wake, and in so doing reduces the number of particles produced opposite to the jet direction. Fig. 10 does not display this well, because in that figure in some panels we have integrated over φ and in others we have integrated over y. To gain a better understanding of this phenomenon, in Fig. 13 we show the energy density profile in the (φ, y) plane, defined as d∆E dφ dy ≡ dp T p 2 T + m 2 cosh y d∆N dp T dφ dy .
For convenience, in the same Figure we show contours corresponding to cones centered on the direction of the high energy parton with fixed cone sizes R. When the cone size is small R < 1, only positive energy from the wake is recovered within the cone. When R increases further and reaches about R ∼ 2, the cone starts to touch the region where the effect of the wake is to deplete the background, represented by a negative energy density change. As a consequence, the accumulated energy recovered by the cone starts to decrease as R increases further in this regime, consistent with what we have seen in the solid curves of Fig. 12.
The details of the shape of the energy ∆E(R) recovered as a function of the cone opening angle R, and in particular the differences between the three solid curves in Fig. 12, reflect details of and differences between the energy density distribution in the three different panels of Fig. 13, one for each of the three calculations. Both for the linearized hydrodynamic model with transverse flow added and for the hybrid model, the positive energy density distribution in Fig. 13 is narrow in rapidity y, and it is recovered first as R increases. The negative energy density contribution is only collected at larger R 2. This is why both the black and magenta solid curves in Fig. 12 drop after R ∼ 2 without further increase.
In contrast, the linearized hydrodynamic calculation leads to a wider y distribution, which implies that not all the positive energy contribution from the wake has been recovered at the angular distance when the effect of the background depletion starts to affect the energy balance. Therefore, the recovered energy in the linearized hydrodynamic calculation without the transverse flow effect added increases first as increasing the cone angle R catches more and more positive energy, then drops with increasing R as negative energy is included, and finally increases again as R increases further, as more positive energy is included. The differences in the role of semi-hard particles in the energy recovered described in the last Subsection prompt us to study the spectrum of particles coming from the backreaction of the medium in further detail. In a realistic simulation, depending on the jet algorithm employed, many of these particles will be a part of what an experimentalist reconstructs as a jet, which will affect the measured fragmentation functions at very low momentum.
Inspired by the definition of fragmentation functions, here we will study the number density of particles from the wake produced within a cone with opening angle R with a given fraction z of the longitudinal momentum of the jet along the jet direction. But, what shall we take as the longitudinal momentum of the jet after it has been quenched? As we have done throughout, for us our "jet" is actually a single high energy parton, not a jet shower. We will assume, again as throughout, that our high energy parton has an initial energy E in = 100 GeV. After it passes through the medium, it has lost ∆E tot . However, we take into account the energy coming from the wake recovered within the cone of opening angle R, namely the quantity ∆E(R) that we plotted in Fig. 12, and define the energy of the "reconstructed jet" as E(R) = E in −∆E tot +∆E(R). For this "reconstructed jet", we can define the longitudinal momentum fraction distribution as ∆f (z) ≡ dp T √ φ 2 +y 2 <R dφ dy d∆N dp T dφ dy ∆f (z) serves as a proxy for the contribution to the in-medium fragmentation function originating from the wake in the plasma.
Results from our calculations of ∆f (z) from all three calculations are shown in Fig. 14 for two different values of R. We have plotted ∆f (z) for z > 0.01, corresponding to particles with p T 1 GeV since the initial energy of the jet is 100 GeV. As we discussed in the Introduction, the oversimplified wake in the hybrid model lacks particles with momentum p T > 2 GeV, which corresponds to z 0.02 here. We see in Fig

VI. CONCLUDING REMARKS
In this paper, we study the wake that a passing high energy parton leaves in the droplet of QGP through which it passes by using linearized hydrodynamics on a background of Bjorken flow. We treat the jet energy and momentum loss as external currents to the hydrodynamic equations and expand the perturbation to these equations describing the jet wake to linear order. By solving the linearized hydrodynamic equations numerically, we have studied the evolution of, and the energy-momentum conservation in, the wake. We observed that viscous effects, even for viscosities as small as η 0 /s = 1/4π, have a significant effect on the evolution of the wake. These effects quickly spread out the hydrodynamic perturbation in the fluid and make the fluid backreaction to an energetic probe small in magnitude. This justifies the linearized hydrodynamics approach we have employed in this paper. By studying particle production from the perturbed hydrodynamic fields at freezeout, we conclude that most of the qualitative features of the backreaction of the medium seen in the complete linear analysis of this paper are well-captured by the simplified approach taken in the hybrid strong/weak coupling model.
We have also studied how radial transverse flow affects the spectrum, azimuthal distribution and momentum rapidity distribution of the particles produced when the medium, including the wake in it, hadronizes at freezeout. We have done so not by solving the full hydrodynamic problem, but rather by taking into account the boost that radial flow exerts on the particles coming from the perturbation of the medium at the time of freezeout. We have seen that this boost has a very significant effect on the particles originating from the backreaction of the medium. As a consequence of the radial transverse flow, these particles originating from the wake in the fluid have a harder spectrum and are more tightly correlated with the direction of the high energy parton than was the case in the simplified hybrid model calculation from previous work. This change in the distribution of the particles coming from the backreaction will have important consequences for jet observables. To get an initial qualitative sense of the extent to which this is the case, we have shown that radial flow effects increase the fraction of the energy lost by the initial high energy parton that ends up (after propagating as a hydrodynamic wake and then freezing out into particles) within a cone of opening angle R around the direction of the high energy parton. This is so for cones with the small values of R used in experimentalists' jet reconstruction algorithms and in fact is so for any R 2. We have also seen that radial flow effects significantly increase the fraction of semi-hard particles with p T > 1 GeV coming from the perturbed medium.
The studies we have performed in this paper cannot be compared directly to experimental measurements, since they lack many important ingredients. For example, we have considered a single high energy parton produced at a specified place with a specified direction, rather than looking at an ensemble of different branching parton showers, with their initial position and direction of origin as well as their pattern of branching drawn from appropriate probability distributions in a Monte Carlo calculation, for example as in the hybrid strong/weak coupling model. Our calculation could also be extended by computing the perturbation on top of a more realistic hydrodynamic background. This notwithstanding, contrasting the results that we have obtained from our linearized hydrodynamic analysis including the effects of radial flow with the much more simplified analysis of particle production from the wake incorporated in the hybrid model is encouraging. Although the hybrid model incorporates the complexities of jet showers and a realistic hydrodynamic background, its oversimplified treatment of the wake has consequences for many jet observables that are sensitive to soft particles within jets [29]. In the detailed analyses of Ref. [29], the hybrid model implementation of particle production from the backreaction of the medium yields a spectrum that is too soft and particle production spread out over too wide an angle to explain experimental measurements of several jet observables in heavy ion collisions at the LHC. We have seen in our analysis that incorporating the effects of radial flow on our linearized hydrodynamic calculation changes the spectrum and distribution of particle production in the direction indicated by the experimental measurements as it yields an increase in the number of semihard particles which are better correlated with the jet direction. While a Monte Carlo study is needed to assess whether the effect introduced by radial flow yields results in quantitative agreement with experimental measurements, the magnitude of the effect we have observed points in the right direction and encourages us to implement these findings in a future Monte Carlo analysis of medium response.
As explained earlier, the physical meaning of the first line in Eq. (A6) is the total energy ∆E deposited into the hydrodynamic system between τ 2 and τ 1 . Writing τ 2 = τ and τ 1 = τ + ∆τ and taking the limit ∆τ → 0, we find τ dx dy dη s J t (τ, x, y, η s ) = d dτ τ dx dy dη s T tτ (τ, x, y, η s ) = dE dτ , which is Eq. (39). Similar equations for the momentum deposition can be worked out and lead to Eq. (40).
We want to emphasize that the condition that the stress-energy tensor vanishes sufficiently quickly at large x, y and η s is important via an example. In the Bjorken flow, the hydrodynamics is boost-invariant, which means the stress-energy tensor is independent of η s In a nutshell, energy is conserved in Bjorken flow. Momentum conservation can be worked out similarly.