A semi-holographic model for heavy-ion collisions

We develop a semi-holographic model for the out-of-equilibrium dynamics during the partonic stages of an ultrarelativistic heavy-ion collision. The model combines a weakly-coupled hard sector, involving gluon modes with energy and momenta of the order of the saturation momentum and relatively large occupation numbers, with a strongly-coupled soft sector, which physically represents the soft gluons radiated by the hard partons. The hard sector is described by perturbative QCD, more precisely, by its semi-classical approximation (the classical Yang-Mills equations) which becomes appropriate when the occupation numbers are large. The soft sector is described by a marginally deformed conformal field theory, which in turn admits a holographic description in terms of classical Einstein's equations in $AdS_5$ with a minimally coupled massless `dilaton'. The model involve two free parameters which characterize the gauge-invariant couplings between the hard and soft sectors. Via these couplings, the hard modes provide dynamical sources for the gravitational equations at the boundary of $AdS_5$ and feel the feedback of the latter as additional soft sources in the classical Yang-Mills equations. Importantly, the initial conditions for this coupled dynamics are fully determined by the hard sector alone, i.e. by perturbative QCD, and are conveniently given by the color glass condensate (CGC) effective theory. We also develop a new semi-holographic picture of jets in the QGP by attaching a non-Abelian charge to the endpoint of the trailing string in $AdS_5$ representing a heavy quark. This leads to modified Nambu-Goto equations for the string which govern the (collisional and radiative) energy loss by the heavy quark towards both hard and soft modes.


Introduction
The experimental heavy ion programs at RHIC and the LHC have brought compelling evidence that the dense partonic form of matter known as the 'quark-gluon plasma' (QGP) created in the intermediate stages of a ultrarelativistic nucleus-nucleus collision behaves in many respects like a 'perfect fluid' -a many-body system with strong interactions. Perhaps the most suggestive data in that sense are those referring to collective phenomena like elliptic flow, which, via their theoretical description in terms of hydrodynamics, point out towards a short thermalization time and a remarkably small value for the ratio between viscosity and entropy density [1,2] (and Refs. therein), as expected at strong coupling. This behavior can be attributed to either large gluon occupation numbers leading to strong non-linear phenomena even when the coupling is weak, or to a genuinely non-perturbative regime, associated with relatively soft gluons, for which the QCD running coupling is indeed quite strong. Most likely, both features combine to yield the observed 'perfect liquid' behavior, and their interplay is expected to vary with time.
The gluonic matter liberated in the early stages of the collision is characterized by large occupation numbers, of order 1/α s , and by a relatively weak value of the QCD running coupling α s = g 2 /4π -precisely because the density is so high. Via the gluon non-linear dynamics, this high density introduces a semi-hard scale, the 'saturation momentum' Q s , which is the typical momentum of the liberated gluons and also the relevant scale for evaluating the running coupling: α s (Q s ) 1. But although the coupling is weak, the mutual interactions of this dense gluonic matter -also called the 'color glass condensate' (CGC) [3,4] are effectively of order one: the factors of α s from the interaction vertices are compensated by the large occupation number ∼ 1/α s of the interacting modes. With increasing time, this matter rapidly expands and thus gets diluted, while at the same time abundantly emitting relatively soft gluons, with energies and momenta well below Q s . These soft gluons have occupation numbers of order one or less, but a relatively stronger coupling, perhaps as strong as α s ∼ 1, and their mutual interactions are expected to represent the driving force towards thermalization. In particular, the semi-hard partons eventually thermalize because they lose energy towards the thermal bath made by the soft gluons. This out-of-equilibrium partonic matter, sometimes referred to as 'glasma' [5], can also include very hard partons, with energies and/or masses much larger than Q s , which are created via binary collisions between the constituents of the incoming nuclei. Such 'hard probes' couple to both the semi-hard and the soft constituents of the surrounding medium, possibly in different ways, and thus can lose energy, or get deflected.
It would be highly desirable to have a unique theoretical description for all the aspects of this complicated non-equilibrium dynamics, but this seems difficult to achieve from first principles. So far, most studies have used a one-sided approach to the problem, which privileges either a weak-coupling scenario, or a strong-coupling one.
In QCD at weak coupling, one has been able to deal with the highly-occupied gluon modes via semiclassical techniques, like numerical solutions to the classical Yang-Mills equations supplemented with fluctuations in the initial conditions (and possibly coupled to non-Abelian Vlasov equations for the harder modes) [6][7][8][9][10][11]. An alternative strategy, which becomes appropriate when the occupation numbers drop down to values much smaller than 1/α s , involves analytic and numerical studies of appropriate transport equations, which include gluon radiation [12][13][14][15][16][17]. Whereas the conventional wisdom, based on the experience with kinetic theory, is that the thermalization time at weak coupling should be parametrically large [12,15,16] -due to the smallness of the cross-sections for partonic interactions (and despite the presence of plasma instabilities [8]) -, more recent numerical studies [11,17] suggest that rapid thermalization is not necessarily in contradiction with the weak coupling picture. For instance, Ref. [11] observed a relatively fast isotropization of the pressure tensor due to the non-linear dynamics of strong classical Yang-Mills fields.
At strong coupling, on the other hand, it is customary to rely on the AdS/CFT correspondence, where thermalization in the conformal N = 4 supersymmetric Yang-Mills (SYM) theory is described as the formation of a black-brane horizon in the radial direction of the dual AdS 5 space-time. This process has been studied via both analytic [18][19][20][21][22][23][24] and numerical [25][26][27][28][29][30][31][32][33][34] methods, but in a purely holographic set-up, that is, starting with initial conditions which describe some non-equilibrium state at strong coupling, whose relevance for the heavy-ion problem is difficult to assert. Perhaps the closest to the heavy-ion scenario is the scattering between two gravitational shock-waves, but even in that case the correspondence to the actual physical situation (the scattering between two dense systems of gluons which are weakly-coupled) is far from being clear, and not truly under control.
The need for combining insights and also techniques from both weak and strongcoupling approaches, with the aim of developing a unified approach, has since long been recognized (see e.g. the discussion in [35]). However, so far it has only been implemented [36] in the form of a hybrid model for 'jet quenching' -the energy loss by an energetic parton propagating through a strongly-coupled plasma. In Ref. [36], the fragmentation of the incoming parton has been described within perturbative QCD, while the interactions between the ensuing parton cascade and the surrounding medium have been modeled by using results from the AdS/CFT correspondence, that is, by assuming the plasma to be exclusively made with strongly interacting modes.
In this paper, we would like to propose a more general strategy for simultaneously describing weak-coupling aspects and strongly-coupling ones, that we shall refer to as semi-holography. This aims at keeping both the (semi)hard and soft degrees of freedom in the medium, together with their various types of interactions: weak interactions among the hard partons, as described by perturbative QCD (including the resummations needed to deal with large gluon occupation numbers), strong interactions among the soft modes, as encoded into a holographic model to be shortly explained, and 'hard-soft' interactions which couple the two types of degrees of freedom in a minimal way. As we shall see, our model involves only a small number of free parameters (two parameters related to hard-soft couplings, plus an additional parameter associated with the hard probes), and therefore can be falsified.
Remarkably, our description can be applied to all the stages of an ultrarelativistic nucleus-nucleus collision, from the initial conditions (which are exactly the same as in the standard CGC approach in pQCD) until hadronization. Albeit it reduces to classical equations of motion, the dynamics generated by our model is quite complex. Due to the specificity of the holographic initial and boundary value problem, the coupled equations of motion for the hard and soft sectors should be solved in a self-consistent manner, via numerical iterations. Our description naturally ends when the effective temperature drops below the confinement scale, meaning at the onset of hadronization.
The holographic part of our model is a minimal version of classical gravity in AdS 5 , which includes the graviton and a massless scalar field (the 'dilaton'). This gravitational theory in AdS 5 may be thought of as dual to an 'infrared conformal field theory' (IR-CFT), which is strongly coupled and provide a model for the soft, non-perturbative degrees of freedom in the QGP phase of QCD. This is completed by the 'hard-soft' interactions, which represent marginal deformations of the strongly-coupled IR-CFT, and in the dual description, provide boundary sources for the graviton and dilaton fields in the bulk.
To summarize, the semi-holographic model contains classical Yang-Mills fields representing the hard sector, ('hard gluons', or 'glasma fields') coupled to a few operators of a strongly coupled IR-CFT, whose spacetime evolution is self-consistently given by the equations of classical gravity in AdS 5 with dynamical sources on the boundary. The couplings between the glasma fields and the IR-CFT operators involve gauge-invariant, local operators of lowest possible scaling dimensions selected according to their quantum numbers. We need only two such operators both in the hard and soft sectors -namely the energy momentum tensor and the 'glueball' operator (the Yang-Mills Lagrangian density). We will derive equations for the coupled evolution of the hard and soft sectors, with initial conditions entirely fixed by the CGC description of the hard modes, and devise a numerical strategy for solving these coupled equations.
Additionally we will also provide a semi-holographic model for the energy loss by a heavy quark. This quark simultaneously couple to the perturbative, glasma, fields -via its non-Abelian, 'color', charge, whose dynamics is modeled by the Wong equation -and to the soft, strongly-coupled, sector -via its virtual cloud of soft quanta, as captured by a trailing string in AdS 5 attached to the non-Abelian charge at the boundary. The Nambu-Goto dynamics of this string in the evolving bulk geometry leads an interesting generalization of the non-Abelian Lorentz force acting on the heavy quark, which encompasses both the respective force due to the classical glasma fields, and the flow of energy and momentum towards the soft degrees of freedom. Using the method in [36], one can extend our model to include the perturbative fragmentation of the heavy quark into quarks and gluons, and thus provide a semi-holographic picture for the evolution of a jet as a whole.
A semi-holographic approach has been previously applied to condensed matter systems, with a motivation similar to ours -namely, to build an effective theory for a large class of emergent non-Fermi liquids close to quantum criticality, integrating both weak and strong-coupling aspects [37,38]. In particular a generalization of Landau's theory for these non-Fermi liquids has been proposed [38], where the entire low energy phenomenology is described in terms of generalized Landau parameters, which may vary from one system to another. However, there are some noticeable differences between such previous approaches and our present one: in the context of condensed matter, the soft degrees of freedom are parametrically denser (by a factor N 2 ) than the hard degrees of freedom in the large-N limit 1 . Accordingly, in that case, one can ignore the back-reaction of the hard degrees of freedom on the soft sector, which greatly simplifies the analysis of the latter.
By contrast, in our semi-holographic model for the QGP, the densities of both hard and soft modes should scale in the same way with N , since they both count the number of gluon-like degrees of freedom 2 (with however have different energies and intrinsic gauge couplings). This feature not only complicates the dynamics in our model, but also introduces some ambiguity in its formal structure: unlike in the corresponding discussion of non-Fermi liquids [38], here it is not possible anymore to distinguish between relevant and irrelevant operators in the large N limit. In view of that, we have opted for a 'minimalist' model involving the smallest number of hard-soft couplings, which is two. One additional argument behind this choice, is the fact that couplings involving "multi-trace IR-CFT operators" should be dynamically suppressed by powers of the hard saturation scale Q s . In Section 3, we will also present a renormalization group argument which gives more support to our model. The plan of the paper is as follows. In Section 2, we will present our semi-holographic model, meaning the full set of coupled equations for the hard and soft degrees of freedom. We will then discuss the initial conditions (for both hard and soft sectors) and explain that they can be unambiguously inferred from the CGC description of the colliding nuclei. We will also describe the iterative numerical strategy for solving the equations. In Section 3, we will present a renormalization group flow argument outlining the set of assumptions under which our model will turn out to be an effective theory for the QGP phase. In Section 4, we will present our semi-holographic model for jets. In the concluding Section, we will mention 1 The parameter N 2 , with N playing the role of the rank of the non-Abelian gauge group in the emerging IR-CFT theory, is related to the inverse of the gravitational constant in the holographic description. 2 In the present, heavy ion, context, this parameter N should be related to the rank Nc of the 'color' gauge group in QCD; see the discussion in the next section.
some key directions which can be readily pursued and some further possible simplifications of our model that may be easier to implement.

The model
In the early stages of heavy ion collisions in QCD, the hard scale is set by the saturation scale Q s of the incoming nuclei. This is so because the partons (mostly gluons) released just after the collision, at a proper time τ ∼ 1/Q s , have both transverse and longitudinal momenta of the order of Q s . The QCD coupling at this scale is generally assumed to be reasonably weak, α(Q s ) 1, for perturbation theory to apply 3 . At times τ 1/Q s , the gluon modes have large occupation numbers, of order 1/α(Q s ), so the non-linear effects are strong even if the coupling is weak. This non-linear dynamics is well accounted for by the classical Yang-Mills equations supplemented with a renormalization group equation for the distribution of their 'sources' (the color charge densities in the incoming nuclei). This picture lies at the basis of a semi-classical effective theory for the initial and early stages of a heavy ion collision, the 'color glass condensate' [3,4], which is derived from pQCD.
However, the classical Yang-Mills equations cannot give an accurate description of the dynamics at later stages, for at least two reasons: on one hand, the partonic system becomes more and more dilute, and hence under-occupied, due to rapid longitudinal expansion; on the other hand, soft gluons are abundantly radiated and they are expected to interact rather strongly with each other since the strength of their mutual coupling is relatively large. These strong interactions should bring the soft modes close to thermal equilibrium, albeit at a rather low temperature to start with, because of their comparatively low energy density. Subsequently, this soft thermal bath will continuously absorb energy from the reservoir made by the hard particles and thus increase its temperature. At the same time, the hard partons will lose energy towards the soft modes, until they reach equilibrium with the latter.
As explained in the Introduction, we here adopt the point of view that the soft gluon modes are strongly coupled among themselves and that their internal dynamics can be effectively described by an infrared conformal field theory (IR-CFT), which admits a holographic dual. More precisely, this effective description applies to the intermediate-energy sector between the saturation (Q s ) and the confinement (Λ) scales, that we assume to be well separated from each other: Q s Λ. We shall require this emergent IR-CFT to be a proper conformal field theory, which exhibits scaling symmetry and Lorentz invariance, and for which one can take the limit of a strong 't Hooft coupling. The Lorentz invariance is of course physically motivatedthis enables us e.g. to generate solutions which describe a boost invariant Bjorken flow for arbitrary initial conditions. The assumption of scaling symmetry, on the other hand, and also the strong-coupling limit, are rather crude idealizations, which are intended to simplify the problem. They amount to neglecting the running of the QCD coupling within the intermediate energy range alluded to above, albeit this running was of course essential to justify the existence of a strong coupling in the first place. Under these assumptions, we can consistently conjecture that the holographic dual of the IR-CFT is a supergravity theory living in an asymptotically AdS 5 space-time.
In the above, we have implicitly assumed that the color gauge group SU(N c ) is the same for the IR-CFT as in QCD. However what we are really assuming here is that the N c of the gauge group of the IR-CFT is parametrically as large as 3 -the N c of QCD, and of course 3 is large enough to a good approximation to enable us to take a large N c limit. There are well known examples in supersymmetric QCD where a new gauge group emerges in the infrared [39] which may or may not be strongly coupled. In our case, the holographic IR-CFT is strongly coupled and is assumed to emerge at intermediate scales in the quark-gluon plasma phase, where the effective coupling is expected to run very slowly. For consistency, we shall use the large N c limit also in the treatment of the 'hard' gluon modes also, which are weakly coupled, and for which we shall keep the semi-classical QCD description. However the 't Hooft coupling λ ≡ g 2 N c for the hard sector is small: λ 1.
We are now in a position to describe our semi-holographic model in more detail. This model combines 'hard' gluon modes, which are weakly coupled but highly populated, and are represented by classical Yang-Mills fields (also known in this context as 'glasma fields'), together with 'soft' and strongly coupled degrees of freedom, which are represented by gauge-invariant operators in the IR-CFT, and whose dynamics can be studied with the help of the supergravity dual. The 'internal' ('t Hooft) couplings in these two sectors are treated as independent quantities, which are both fixed : in the hard sector, the coupling is small, λ 1 as mentioned above, whereas in the IR-CFT, the respective coupling is infinitely strong and disappears from the problem in the classical gravity approximation. These two types of degrees of freedom must be coupled with each other via operators which become irrelevant in the 'hard' limit Q s → ∞ -as in this limit the dynamics of the hard modes should become insensitive to the soft sector. Furthermore, the full hard and soft sectors must have a consistent large N c limit.
Since the relevant operators of the IR-CFT are gauge invariant, the hard modes can couple to them only via gauge-invariant composite operators. A minimal model involves only two marginal operators from the IR-CFT, the energy-momentum tensor of the soft modes that we denote as T µν and a scalar operator, denoted as H, which can be thought off as the glueball operator in the IR-CFT. The classical gravity fields dual to these operators are the graviton and a massless dilaton respectively. They obey the standard equations of classical gravity in the dual AdS 5 geometry, to be presented below. The 'partonic' operators which couple to these IR-CFT operators can be inferred from symmetry arguments, together with the requirement that the hard-soft couplings be suppressed by inverse powers of Q s . As the expectation values of the soft sector operators will be vanishing at initial times (the soft sector can be produced only later by radiation from the hard sector), they can be expected to remain small in units of Q s throughout the entire evolution, due to dilution by the longitudinal expansion. Therefore the coupling to other gauge-invariant operators of the soft-sector of higher scaling dimensions will be suppressed by powers of Q s .

The coupled equations of motion
The minimalist approach to coupling the IR-CFT energy-momentum tensor T µν to the weakly interacting glasma fields is a tensorial coupling 4 to t µν , the energy-momentum tensor of the Yang-Mills fields. The latter has the standard expression in matrix notation appropriate for the adjoint representation, where the factor N c has been introduced by the trace normalization in the adjoint representation -Tr(T a T b ) = N c δ ab , with a, b = 1 · · · N 2 c − 1 representing the color indices for SU(N c ) gauge group. We use the standard convention for the (gauge) covariant derivative: where g is the effective Yang-Mills coupling for the glasma fields as discussed above.
Furthermore, the scalar operator H of the IR-CFT will naturally couple to the 'glueball' operator for the glasma fields, i.e. the density of the Yang-Mills action: Based on these considerations and dimensional arguments, we are led to the following minimalist action describing the dynamics of the weakly-interacting glasma fields in the presence of the strongly-interacting soft radiation: Above, β and γ are two dimensionless parameters in this effective theory, which at large N c must scale like 1/N 2 c . In this case, all the terms in the above action scale in the same way at large N c , since the two IR-CFT operators H and T µν are expected to scale like N 2 c in the large N c limit. On the other hand, the expectation values of these IR-CFT operators are also expected to remain much smaller than Q 4 s for reasons discussed above, hence the 'hard-soft' couplings in Eq. (2.3) should act as a small perturbation on the dynamics of the hard modes. Also, with a slight abuse of notations, we have used the notations H and T µν not only for the respective IR-CFT operators, but for the corresponding expectation values in the IR-CFT state to be determined self-consistently. We will continue to use these notations, it should be clear from the context whether we are referring to the IR-CFT operators or their expectation values.
The coupling to the IR-CFT modifies the equations for the glasma fields (inside the forward light-cone whose tip is at the event of collision; see also the discussion in Sect. 2.2 for more details) which can be readily obtained from Eq. (2.3) and read As announced, the dynamics of the strongly coupled IR-CFT will be studied using holographic methods. Via the 'hard-soft' couplings in Eq. (2.3), the glasma degrees of freedom act as sources for the IR-CFT, which is thus marginally deformed. In the dual, holographic, description, this implies that the dual fields in the bulk acquire non-trivial boundary values, so that the non-normalizable modes of these dual fields are non-vanishing. These boundary values must equal to (β/Q 4 s )h for the bulk dilaton field (dual to the IR-CFT operator H) and respectively (γ/Q 4 s )t µν for the bulk graviton (dual to the IR-CFT operator T µν ).
Let us denote the massless dilaton field in the bulk as Φ and the five-dimensional metric as G M N . We use the Fefferman-Graham coordinates for the bulk spacetime, where x M = (z, x µ ), with z denoting the holographic radial coordinate and x µ the four fieldtheory coordinates. In these coordinates, the boundary of AdS 5 is located at z = 0. Then the near-boundary behavior of the bulk fields Φ and G M N reads where H and T µν denote the respective expectation values in the dual IR-CFT state. The additional tensor P µν occurring in the coefficient of the z 4 term in the above expansion of G µν is a well-defined functional of the leading term g (0)µν ≡ η µν + (γ/Q 4 s )t µν , which appears because this leading term itself (which, we recall, represents the metric in the dual IR-CFT theory) is different from the Minkowski metric η µν and hence it describes a curved space-time. (We shall shortly return to a discussion of the physical meaning of this nontrivial boundary metric g (0)µν .) This tensor P µν can be constructed following the method described in Ref. [40], which yields where g (2)µν is the coefficient of z 2 in the expansion of G µν , namely In the above expressions, all raising and lowering of indices have been defined with the metric g (0)µν or its inverse. Since, however, g (0)µν = η µν + (γ/Q 4 s )t µν , one sees that P µν can be also viewed as a functional of t µν in the flat, Minkowski space-time.
Thus by expanding the self-consistently determined solution of the five-dimensional classical gravity equations, to be specified soon, in a power series of the radial coordinate near the boundary, one can obtain both H and T µν . In Eq. (2.5), l is the 'AdS radius' and determines the bulk five-dimensional cosmological constant via Λ 5 = −6/l 2 , and G 5 denotes the five-dimensional Newton's constant. Note that l and G 5 do not appear separately in any physical quantity 5 , but only in the dimensionless combination 4πG 5 /l 3 , which according to the AdS/CFT dictionary should be identified with 2π 2 /N 2 c . Recalling that both β and γ scale like 1/N 2 c , we conclude that all the terms shown in the expansion in Eq. (2.5) have the same scaling in the large N c limit: they are all of O(1).
The five dimensional classical gravity theory giving equations for the dynamics of the massless bulk fields Φ and G M N can be taken to be those of Einstein's gravity with a minimally coupled massless dilaton and a negative cosmological constant, because a twoderivative action in the bulk is expected in the strong coupling limit. By also including a properly chosen dilaton potential V (Φ), as e.g. in [41], one could mimic the running of the QCD coupling in the 'soft' sector. For simplicity, here we shall limit ourselves to the case of a fixed coupling, so there is no potential for the dilaton. Then the bulk equations of motion take the standard form These equations must be solved with the boundary conditions as z → 0 that can be read from Eq. (2.5) together with initial conditions appropriate for the problem of heavy-ion collisions, to be described in the next section. These constraints together with the requirement that there must be a regular horizon in the bulk ensure the unicity of the solutions.
In the standard holographic dictionary, the boundary metric, which is the leading term in the near-boundary expansion of G µν in Eq. (2.5), that is, is the metric in which the degrees of freedom of the boundary gauge theory 'live'. In our case, this 'boundary gauge theory' is of course the IR-CFT, which is not a fundamental theory by itself, but an effective theory for the strongly-coupled, soft, modes of the glasma. The fact that this effective theory 'lives' in a non-trivial metric should not be meant to imply that these non-perturbative degrees of freedom truly live in some curved space-time.
Rather, it simply expresses the fact that the IR-CFT is 'deformed' by its coupling to the hard sector, meaning that hard and soft d.o.f. exchange energy with each other. As we now explain, this interpretation naturally emerges from the Einstein equations in the bulk. Specifically, in the limit z → 0, one of the constraints in the Einstein equations imply the following conservation law for energy and momentum: where the space-time indices are raised/lowered with the help of the boundary metric (2.9) or its inverse, e.g., 11) and the covariant derivatives are written with respect to this metric as well. Recalling that the Christoffel symbols for the Levi-Civita connection Γ associated to this metric are built with the energy-momentum tensor of the hard modes, (in the r.h.s. of this equation, all lowering and raising of indices are done with respect to the flat metric η µν and its inverse), one sees that the conservation law in Eq. (2.10) can be rewritten in terms of ordinary derivatives in flat, Minkowski, space-time and can be properly interpreted as the standard conservation of energy and momentum, but in the presence of an external 'driving force' representing the hard degrees of freedom; that is, This is precisely the expected conservation law for a generic field theory (here, our IR-CFT) which is 'deformed' by the irrelevant operators shown in Eq. (2.3). Independently of this holographic context, this law can also be obtained as a Ward identity via standard manipulations on the path integral for the IR-CFT. Notice that, although exact (in an operator sense), this equation (2.13) is not sufficient to determine the dynamics of the energy-momentum tensor in the IR-CFT (for given hard 'sources' h and t µν ), because there are only 4 such equations constraining 10 independent degrees of freedom. As a matter of fact, there is still another constraint -an equation for the trace of T µν as computed in the boundary metric, which should be reinterpreted too as an equation in flat space involving hard 'sources' -, but the number of such constraints is still too low to fully fix the dynamics, as expected. One needs to solve the complete bulk dynamics, meaning (2.8), to obtain the evolution of T µν and H. The previous discussion also points towards a rather subtle point, that we shall now try to elucidate: since, as we have just seen, the IR-CFT sector formally lives in a curved space-time with metric g (0)µν (x), whereas the Yang-Mills sector lives in the physical spacetime with Minkowski metric η µν , there is a priori an ambiguity in the construction of the tensorial coupling between the energy-momentum tensors of these two sectors. In writing the action (2.3), we have chosen this coupling in the form t µν T µν , which is of course 6 the same as t µνT µν , withT µν ≡ η µα η νβ T αβ , but is not the same as t µν T µν , with T µν the covariant tensor which enters the near-boundary expansion of G µν shown in (2.5). (The relation between the covariant and contravariant components of this tensor is rather given by Eq. (2.11).) Our above choice to minimally couple to contravariant tensor T µν to the corresponding tensor of the hard modes -that is, to define the hard-soft coupling as t µν T µν , rather than t µν T µν -is truly a prescription, which should be viewed as a part of our model and is ultimately guided by physical considerations: from the viewpoint of the IR-CFT theory, the 'source' t µν plays the same role as a deformation δg µν in the metric, and hence it is most naturally considered with covariant indices.
To conclude, the glasma fields F a µν and the expectation values H and T µν of the IR-CFT operators must be simultaneously determined by self-consistently solving the generalized Yang-Mills equations (2.4) and the classical gravity equations (2.8), with the boundary conditions (2.5) for the latter. For the corresponding solutions to be uniquely determined, we still need to specify the initial conditions at τ = 0 and the regularity condition in the bulk. These additional ingredients will be discussed in the next two subsections, together with more details on the numerical implementation of the self-consistency procedure.
Notice that the time scale for thermalization and the onset of a hydrodynamic description will be affected by both the glasma and the holographic degrees of freedom. So, our scenario cannot be expected to predict universal values for the thermalization time and the transport coefficients, as is generally the case for the purely holographic models at strong coupling. (In particular, the non-linear dynamics of the glasma fields can by itself lead to fast isotropization [11].) All observables will be determined in terms of the two free parameters of the model, namely the couplings β and γ in the action (2.3), together with the dimensionless quantity 4πG 5 /l 3 = 2π 2 /N 2 c of the classical gravity theory in the bulk.

Initial conditions and the iterative algorithm
A main virtue of our present semi-holographic model is that it provides an unambiguous initial-condition problem (including for the holographic sector), which is entirely determined by the pQCD picture for the colliding nuclei -more precisely, by the associated CGC-description in terms of color charge distributions [3,4] (see also below). This should be contrasted to purely holographic approaches, which requires an ad-hoc initial condition for the gravitational problem, that cannot be directly related to properties of the colliding nuclei.
It is first of all clear that the initial conditions in the hard sector, i.e. for the classical glasma fields, will not be affected by the holographic degrees of freedom. Using causality and the uncertainty principle, one can argue that the expectation values of T µν and H (the soft sector operators) must vanish at proper time τ = 0 (the time where the collision is initiated) : indeed, it takes some time τ ∼ 1/p ⊥ 1/Q s to emit soft gluons with transverse momenta p ⊥ Q s . Hence, when τ → 0, the modified glasma equations (2.4) reduce to the standard Yang-Mills equations and the initial conditions for our glasma fields are the same as in the standard CGC approach at weak coupling. These initial conditions have been worked out in [42,43] and will be briefly discussed below.
Before revising the standard CGC initial conditions for the glasma fields, let us also consider the corresponding problem on the gravitational side. As is well known, the gravitational problem in an asymptotically AdS 5 spacetime cannot be reduced to an initial-value problem alone. Rather, one also need to specify the boundary (z → 0) metric and the boundary value of the dilaton field at all times -to ensure a well-defined evolution that leads to a space-time with a smooth future horizon. In our case, these boundary data are specified by (2.5) in terms of the glasma fields. However, in order to know the latter at all times, we need to solve the modified glasma equations (2.4), which in turn involve the 'soft' operators T µν and H, as determined by the solutions to the gravitational problem. Hence, we are facing a self-consistent problem, whose solutions calls for an iterative procedure. In what follows, we will lay out a well-defined procedure in that sense.
Let us first review the standard CGC initial conditions for the glasma fields [42,43]. To that aim, we recall the standard glasma equations, which do not include the soft degrees of freedom, but where the color sources representing the colliding nuclei are manifest 7 : (2.14) We are using light-cone coordinates: is the color current of the first nucleus, which is a 'right mover' (it propagates towards increasing x 3 at nearly the speed of light). The corresponding color charge density ρ a (1) (x − , x ⊥ ) is independent of x + by Lorentz time dilation and is localized near . A similar discussion, with x + ↔ x − , applies to the second nucleus, which is a left mover and has a color current J ν a (2) The color charge densities are randomly distributed, independently for each nucleus, according to a gauge-invariant probability distribution (the 'CGC weight function'), which in the simplest approximation (known as the 'McLerran-Venugopalan model') is taken to be a Gaussian: In this equation, the indices m, n = 1, 2 refer to the two nuclei and the 2-dimensional density of color charge squared g 2 µ 2 is proportional to the saturation scale Q 2 s [44][45][46]. Eq. (2.14) is written in the Schwinger gauge 2) is not rotated by the classical field 8 . The solution to this equation outside the forward light-cone -meaning prior to the collision (t < 0) and in the region which is causally disconnected from the collision at t > 0 -is the superposition of the individual fields of the two nuclei, which outside the light cones are 2-dimensional pure gauges : where the upper index i = 1, 2 denotes the transverse component and with the P symbol indicating the path ordering of the color matrices ρ (1,2) = ρ a (1,2) T a in the exponent w.r.t. the relevant longitudinal coordinate, x − or x + .
After the collision, the solution to Eq. (2.14) in the future light cone (t > |x 3 |) is boost-invariant, meaning that it depends upon the proper-time τ , but not on the rapidity η. These variables τ and η are defined as (for t > |x 3 |) More precisely, the solution (in the gauge , with the functions A η and A i obeying the following initial conditions at τ = 0 : where the square brackets refer to the commutator of color matrices. As previously mentioned, in our semi-holographic approach, the initial conditions (2.19) also applies to the complete equation of motion for the glasma fields, i.e. Eq. (2.4) which includes the back reaction of the soft modes. For this equation to be fully specified though, we also need the expectation values H and T µν of the IR-CFT operators, which in turn are determined by solving the gravitational equations with boundary conditions specified by the glasma fields, cf. Eq. (2.5). As discussed earlier, this gravitational problem is not just an initial-value problem: it requires specifying the boundary data at all times. To cope with this, we propose a self-consistent procedure in which the back reaction of the soft sector on the glasma fields is formally treated as a 'small perturbation', which is taken into account via successive iterations -to be repeated until one achieves convergence.
Specifically, in the first iteration, one solves the glasma equation without the holographic d.o.f. in the forward light-cone, that is the homogeneous version of Eq. (2.14) with the initial conditions (2.19). This numerical solution, which is well known in the literature [47,48], is then used to compute the boundary conditions for the classical gravity fields, at all times. This requires performing the average over the distribution (2.15) of the color charge densities in the colliding nuclei, which in practice is done by randomly selecting the initial color charge densities according to the Gaussian distribution (2.15). Strictly speaking, the quantities t µν and h thus obtained exhibit logarithmic divergences as τ → 0 and require an ultraviolet cutoff. However, as argued in Ref. [49], this divergence has no influence on the soft sector because t µν and h become independent of the UV cutoff at later times τ 1/Q s , which is when we actually need them (since, as previously explained, it takes some finite time τ 1/Q s to emit the soft fields). To avoid such complications with the limit τ → 0 in the numerical simulations, we propose to smoothly switch on the interactions between the hard and soft sectors, by replacing the previous couplings β and γ in Eq. (2.3) with their time-dependent versions where α 0 , β 0 and γ 0 are new, dimensionless parameters. This in particular implies that the boundary sources for the gravitational problem will acquire a similar time-dependence: at the boundary r = 0. Similarly the modified glasma equations (2.4) should now read As clear from (2.21), on the gravitational side, due to vanishing of the hard-soft couplings at the initial time τ = 0 + , the boundary metric is just the flat metric η µν and the boundary value of the dilaton field is vanishing at the initial time. Therefore, we can consistently put the initial conditions on the gravity side to be the vacuum pure AdS 5 solution with a vanishing dilaton. Following [33], we shall employ the ingoing Eddington-Finkelstein coordinates, in which the bulk space-time metric can be written as where X denotes collectively all bulk coordinates, r is the new radial coordinate, τ is the proper time defined in (2.18), and x m denotes the three spatial coordinates, namely, the 'rapidity coordinate' η defined in (2.18) and the transverse coordinates x i with i = 1, 2 9 . In these coordinates, the boundary is at r = ∞. The components of the boundary metric g µν as written in (2.21) are related to the boundary limit of the components of G µν of the new metric via Similarly in the second equation of (2.21), the right hand side refers now to the limit r → ∞ of the dilaton field Φ in the new coordinates. In order to ensure the "vacuum initial conditions" for the soft sector, we will require the following initial conditions at τ = 0 + 10 (2.25) 9 The relation between the radial Fefferman-Graham coordinate z and the radial Eddington-Finkelstein coordinate r is r = l 2 /z + O(z 0 ) near the boundary. The relation between the two time coordinates τFG and τEF is τEF = τFG − z + O(z 2 ), etc. near the boundary. All the field-theory coordinates x µ in both these coordinate systems agree at the boundary but they differ in the interior (see [50] for instance). We do not use different notations for the field-theory coordinates in these two coordinate systems here for simplicity. 10 We recall that the components Tµν specify the subleading term, of order z 4 , in the near-boundary expansion of the metric Gµν in Fefferman-Graham coordinates, cf. Eq. (2.5). The components of Tµν can also be read from a similar asymptotic expansion of the Eddington-Finkelstein coordinates (see [33] for instance).
Notice that we need to specify the spatial components G mn entirely (meaning for all values of r), while for the functions A and F m we need only the 'subleading normalizable modes', which are T τ τ and T τ m respectively. This is because A is an auxiliary variable whose timederivatives do not appear in the classical gravity equations. Furthermore, after separating the classical gravity equations into equations for the dynamical evolution and constraints, it follows that time derivatives of F m also do not appear in the actual time evolution equations (see [33] for more details). The radial equations for A and F m can be integrated with the boundary conditions shown in (2.21) and (2.25) 11 . The constraints of classical gravity equations, namely (2.10), give the time evolution of T τ τ and T τ m .
The Chesler-Yaffe algorithm [33] can then be employed to solve the gravitational equations (2.8) with the initial conditions and boundary conditions specified as above. This algorithm ensures that the solution is smooth at the future horizon thus protecting causality in the dual field theory. The key to obtain the regularity at the future event horizon is to enforce an apparent horizon in the bulk geometry at a fixed value of the radial coordinate, which in our case will be put at r = l 2 Λ QCD at all times. We will discuss more about this in the next subsection. This regularity condition fixes the values of the remaining integration constants, thus fully specifying the solutions to the classical gravity equations.
From the above regular solution of gravity, we can extract T µν (thus T µν via (2.11)) and H by doing the near-boundary expansion of G M N and Φ respectively (recall (2.5)). We can then use it in (2.22) to correct the solution for the glasma fields in the next step of the iteration. The corresponding glasma solution should be used to recalculate the sources at the boundary (2.21) in the gravitational problem, leading to new expectation values for the IR-CFT operators. As announced, this iterative procedure should be repeated until we obtain convergence for all quantities (hard and soft) and at all times.
The iterative procedure is summarized below.
• In the first step of the iteration, we solve the glasma equations (2.22) on the forward light cone with the initial conditions (2.19) and the expectation values of the softsector operators T µν and H set to zero.
• The numerical solution of the glasma equations is then used to evaluate the sources (2.21) in the gravitational problem.
• With the above boundary sources and the initial conditions (2.25), the gravitational equations (2.8) are solved numerically using the Chesler-Yaffe procedure ensuring that the solution is regular at the future horizon.
• The above are then used in the glasma equations of motion (2.22) to correct for the solution of the glasma fields obtained in the first iteration. The same initial conditions (2.19) are used to obtain the new solution for the glasma fields.
• The corrected glasma solution is then used to recalculate the boundary sources (2.21) in the gravitational problem.
• The gravitational equations for the bulk geometry and bulk dilaton are then solved again with the above boundary conditions and the same initial conditions (2.25).
• The corrected time-dependent expectation values of T µν and H are then extracted from the gravitational solution and used as inputs to correct the glasma solution.
• The iteration between the solution for the glasma fields and the solution for the gravitational equations continues until we reach convergence for both. Thus we obtain the final solution for both hard and soft fields.
We emphasize that the initial conditions for the glasma fields (namely (2.19)) and for the bulk metric and the dilaton (namely (2.25)) remain the same at every step in the iteration. Accordingly, the only physical ingredient in fixing the initial value problem is the distribution of the color sources in the colliding nuclei, as given by the MV model (or, more generally, by the CGC effective theory [3,4]). As previously noticed, both the glasma equations and gravitational equations with similar degree of complexity have already been solved in the past, via numerical methods. Therefore, it seems that our algorithmic procedure can be readily implemented with currently available numerical codes. It remains to be seen if, in practice, the convergence of our iterative procedure can be achieved fast enough to allow for explicit studies of thermalization and for phenomenological predictions.

The apparent horizon and hadronization
The Chesler-Yaffe algorithm [33] requires cutting out the gravitational metric behind the apparent horizon, which by exploitation of a residual gauge symmetry is fixed to a constant value of the radial coordinate. This residual gauge symmetry is the freedom of changing the radial coordinate from r to r via r = r + λ(τ, η, x i ), which keeps the form of the metric (2.23) invariant. As λ does not depend on the radial coordinate, this freedom can be always utilized to keep the apparent horizon at a fixed value of the radial coordinate. In practice this makes the variable λ dynamical and also requires us to provide its value at the initial time. The natural initial value for λ in our case is λ = 0 as this makes the spacetime metric (2.23) match to a pure AdS metric on a standard constant time slice, at the initial time.
As the actual event horizon is 'above' the apparent horizon (meaning closer to the boundary), this is sufficient to ensure that the solutions for both the bulk metric and the dilaton are smooth at the event horizon, thus satisfying the required regularity condition. The formation of the event horizon in the bulk signifies entropy production and then thermalization to an expanding perfect-fluid like state dual to an accelerating black hole [18].
In our case, the apparent horizon has a specific physical significance by itself. The radial coordinate acts as the scale of a renormalization group flow (see next Section for more details), thus we can think of the apparent horizon as a (geometric) infrared regulator for the soft modes. Indeed, the QCD modes with momenta below Λ QCD (the confinement scale) are not correctly accounted for by our IR-CFT and its holographic dual. To that aim, we shall demand that in some precise set of coordinates we should get an apparent horizon in the bulk at r = l 2 Λ QCD . This will give a concrete realization of the expectation that degrees of freedom below the confinement scale should play no role in the dynamics leading to thermalization.
The effective temperature at late time can be read off by matching the late time behavior of T µν to a hydrodynamic asymptotic series expansion in the proper time τ [18]. At later times, the fireball will cool down due to expansion, meaning that the effective temperature will drop below Λ QCD . This time will define the time of hadronization. The latter can be modeled to in our approach by changing the apparent horizon r = l 2 Λ QCD into a hard wall cut-off immediately when the effective temperature equals Λ QCD . The hard wall cut-off denotes boundary conditions which render the bulk Hamiltonian Hermitean. This has previously been used as an AdS model of confinement in QCD [51][52][53]. The bulk excitations of both the metric and the dilaton can be decomposed into eigenfunctions each of which represent specific hadronic excitations 12 . This will be legitimate because the glasma fields will have dissipated into soft excitations and will be negligibly small -and thus will play no role in the hadronization process.

A heuristic derivation using renormalization group flow
It is generally admitted that the radial dimension in the gravity set-up somehow captures the scale of a special renormalization group (RG) flow in the dual field theory, however the way how this actually happens is one of the outstanding problems regarding a fundamental understanding of the AdS/CFT correspondence.
One can think of an abstract definition of the strong interaction limit as opening of a large gap in the scaling dimensions of the operators 13 . Furthermore, the operators which survive this limit form an algebra generated by a few, single-trace, operators, namely the energy-momentum tensor T µν , conserved currents J µ , and relevant condensates H. All the other operators in the strong interaction limit are multi-trace operators like T µν T ρσ .
Another pre-requisite for the field theory to have a holographic classical gravity dual, is the presence of a parameter N , such that in the large N limit, the expectation values of 12 To be more realistic we need to also add gauge fields in the bulk which will capture the hadronic excitations with flavour quantum numbers [53]. These will be also needed to understand the role of chemical potentials in the thermalization process. We can ignore these complications because the chemical potentials are expected to be small. 13 This has to be understood in the sense of the parametric dependence of the anomalous dimensions of the various operators upon a field-theory coupling, such as the 't Hooft coupling in Yang-Mills theory. For most operators in the theory, the anomalous dimensions will receive parametrically large quantum corrections, therefore at infinitely strong coupling the corresponding operators decouple. The operators which survive in the strong coupling limit are those which are protected by symmetries and for which the scaling dimensions remain bounded even when the coupling becomes arbitrarily large.
where the expectation value refers to any state of the theory. We note that, despite this factorization property, the single trace operators will mix with multi-trace operators under the RG flow, thus giving rise to non-linearities. This is the fundamental reason why the dual classical gravity theory is non-linear. Let us give here a general argument, to be developed in detail in a further publication [54], about how the classical gravity theory can emerge from the dual field theory via a special renormalization group (RG) flow. We shall first discuss the pure holographic setup and then generalize to the case of semi-holography. Our argument will be restricted to the hydrodynamic sector describing the dynamics on large space-time scales, where the evolution is governed by local conservation laws. Under this assumption, we would like to argue that the gauge/gravity correspondence can be understood as a reconstruction of a special, highly efficient, RG flow, which can only exist in the large N limit (a condition which is necessary but in general not sufficient).
Consider projecting the degrees of freedom in the field theory on a certain slowly varying subsector via a projection operator P(Λ), where Λ labels the typical scale for space-time variation for fields within this subsector. We also consider the following projection of an observable like the energy-momentum tensor into that subsector: The result of such a projection is a scale-dependent, coarse-grained, observable. This observable will not satisfy the usual microscopic Heisenberg equation of motion. For instance, unlike the microscopic energy-momentum tensor which satisfies the standard conservation equation ∂ µ T µν = 0, its coarse-grained version T µν (Λ) will satisfy a more complicated equation reflecting the fact that the hard degrees of freedom that have been projected out will generate driving forces for the slow moving sector. For simplicity, we have assumed that T µν is the only single-trace operator. The non-linear contributions on the right hand side are due to the mixing of T µν (Λ) with the multi-trace operators (T µν T ρσ )(Λ), etc. under the RG flow. The main point for constructing a "highly efficient RG flow" is to construct a sequence of special projection operators P(Λ) such that the above equation takes a very special form, namely This form is special in that it can be recognized as a covariant conservation law in the effective, scale-dependent, metric g µν (Λ), defined as Indeed, Eq. (3.4) can be rewritten as where ∇ (Λ) is built from the metric (3.5). Thus the equation for T µν (Λ) takes exactly the same form as the usual conservation equation ∂ µ T µν = 0, but in a new, effective, metric. Therefore the notion of "highly efficient RG flow" is that under this RG flow the equations obeyed by the projected operators take the same form as in the ultraviolet, but with an effective background metric which varies with the RG scale in such a way to absorb the contributions from multi-trace operators to the respective equations. Besides this effective metric, the equations for the projected operators can also include additional sources coupling to single-trace operators, which are scale-dependent as well 14 .
One can reconstruct the dual classical gravity theory from this "highly efficient RG flow" as follows. After taking the expectation value of Eq. (3.5) (assuming the large-N factorization of multi-trace operators) and identifying the radial coordinate z in the bulk spacetime with the scale Λ via z = Λ −1 , the effective metric g µν (Λ) becomes essentially the induced metric on the hypersurface z = Λ −1 in the bulk spacetime. Here we have assumed that the boundary (the UV of the dual theory) is at z = 0, so that the radial coordinate is the inverse of the scale. Nothing much changes if we do a radial reparametrization to put the boundary at r = ∞, using as for instance r = l 2 /z, except that the scale Λ now should be related to r by r = l 2 Λ. The effective equation (3.6) is then equivalent to the so-called momentum constraints of the equations of classical gravity, evaluated on this hypersurface, with T µν (Λ) being related to the extrinsic curvature of the hypersurface.
In order to reconstruct the full bulk spacetime metric, one needs to decipher the automorphism group of the renormalization group flow which is related to the lift of the ultraviolet Weyl symmetry to an arbitrary scale. It can be shown that the latter completely encodes the information of the corresponding gauge fixing of the diffeomorphism symmetry in the bulk [54]. Using this information, one can recover the complete bulk spacetime metric from the renormalization group flow with an additional boundary condition given by the boundary metric, which corresponds to the actual metric on which the dual field theory lives [58,59]. This leads to the understanding that the holographic correspondence is equivalent to an existence of a special renormalization group flow in the field theory that is well defined by a set of principles.
The above discussion applies for purely holographic field theories. Our semi-holographic model can be motivated via three simple postulates: 14 The idea that the sources coupling to single-trace operators become dynamical when the multi-trace operators are integrated out leading to classical gravity in one higher dimension has been explored before in the path-integral viewpoint before [55][56][57]. The approach described here takes the Heisenberg picture point of view.
• A highly efficient RG flow in the sense described here can be constructed only for a subset of soft degrees of freedom.
• The effective background for these subset of degrees of freedom in the ultraviolet (i.e. initial conditions for the RG flow of this subset of degrees of freedom) are given by expectation values of operators projected to the hard sector.
• An effective intermediate scale (like the saturation scale Q s here) controls the way the hard sector determines the effective background for the soft sector in the ultraviolet.
These postulates motivate our model, more particularly the action (2.3) entirely. A more precise formulation is beyond the scope of this paper. Let us finally notice that we have already used the above perspective in reinterpreting equations like (2.10), which were a priori obtained within holography and express conservation laws in a fictitious, curved, space-time metric (similarly to Eq. (3.6)), as ordinary equations of motion in flat space-time, but with additional sources describing the coupling to the hard sector, cf. Eq. (2.13). As compared to Eq. (3.4), the right hand side of Eq. (2.13) does not explicitly contain terms non-linear in the soft operator T µν , but only its couplings to hard sector observables, like t µν . Yet, such non-linear effects in T µν are eventually generated via the backreaction of the soft sector on the hard fields, as encoded in the self-consistent solutions to the coupled equations of motion for the two types of degrees of freedom.

A semi-holographic model for jets
We propose the following semi-holographic picture for a 'jet' -more precisely, a very heavy quark in the fundamental representation of the QCD gauge group SU(N c ) -propagating through the non-equilibrium quark-gluon plasma. We represent the hard part of the quark wavefunction as a classical pointlike particle, which carries a non-Abelian color charge Q a and interacts perturbatively with the Yang-Mills fields in the glasma. The cloud of virtual soft particles surrounding the hard core can be thought of as a usual probe string trailing in the bulk. The latter follows Nambu-Goto equations, which are generalizations of the geodesic equation for a one dimensional object and capture the influence of the strongly coupled degrees of freedom on the soft sector on the quark wavefunction.
Let us denote the world-sheet coordinates of the string as σ α = (τ, σ). The string embedding in the bulk spacetime is given by X M (τ, σ), where X M = (r, x µ ) are the coordinates of the bulk spacetime. The world-sheet time coordinate τ parametrizes the boundaries which are world-lines. The boundary σ = 0, in particular, will correspond to the end point of the string fixed on the boundary of the bulk AdS 5 spacetime, which is given by r = ∞ (or z = 0 in Fefferman-Graham coordinates). (We assume an infinitely heavy quark, for simplicity). We will use Eddington-Finkelstein coordinates in the bulk for the rest of this Section.
The semi-holographic model for the dynamics of the open string is encoded in the following action: where T 0 is the string tension and h αβ is the induced world-sheet metric, which reads with the string-frame metric G It is not hard to see that T 0 enters the equation of motions only in the combination T 0 l 2 , where l is the radius of the ambient anti-de Sitter space. Therefore, we have one additional dimensionless parameter, namely T 0 l 2 , for describing jet dynamics.
In order for the action to be gauge-invariant, the non-Abelian charge Q a at the end of the string should follow the following condition (below, which is recognized as the Wong equation [60] adapted to the problem at hand (see e.g. [61]). For the physical interpretation and also for practical calculations, it is useful to chose the world-sheet coordinates as τ = t, σ = r, (4.5) where t and r are the time and radial coordinates of the embedding spacetime. In this case, the bulk coordinates X M (τ, σ) become t, r, x i (t, r) . We recall that the radial coordinate r extends from ∞ up to the apparent horizon at l 2 Λ QCD . In this gauge, the world-sheet action becomes The last piece in the action is then recognized as the standard interaction between a (classical) color current and the non-Abelian gauge field, that is, After adding the string action (4.6) to the glasma action (2.3), it is clear that the color current j µ a (x) of the heavy quark counts as an additional source in the r.h.s. of the Yang-Mills equations (2.4). By assumption, the heavy quark represents only a small perturbation of the medium, so the change in the glasma field induced by this current can be computed in the linear response approximation. The small correction δA µ obtained in this way, which is linear in j µ , is important for our present purposes, in that it encodes a part of the energy loss by the heavy quark -namely, that part due to radiation which is induced by its interactions with the hard, glasma fields. In a purely perturbative set-up at weak coupling, a similar calculation provides the medium induced radiation a la BDMPSZ [62][63][64]. Note that the interactions relevant for this mechanism involve both the color precession of the color current j µ a (x) of the heavy quark, as described by Eq. (4.4), and the non-linear effects in the Yang-Mills equations (2.4), which describe the scattering between the radiated gluon (the induced field δA µ ) and the hard modes in the glasma (the background field A µ ).
In the present set-up, there are two additional mechanisms leading to energy loss by the heavy quark. First, there is collisional energy loss associated with the scattering off the hard partons; in the present approximations, this is described as the acceleration of the heavy quark by the Lorentz force generated by the glasma fields. Second, there is another type of medium-induced radiation, which is generated by the interactions between the heavy quark and the soft, strongly-coupled, modes; in our model, this is holographically represented as the energy flow down the string world-sheet. To compute these various effects, one needs the trajectory x i (t) of the heavy quark (the string endpoint at the Minkowski boundary) and, more generally, the string world-sheet X M (τ, σ) in the bulk. For that purpose, we need to solve the Nambu-Goto equations with appropriate boundary conditions.
To construct these equations and the associated boundary conditions, it is preferable to work with manifestly geometric notations, which are reparametrization-invariant. Hence, we return to our generic world-sheet coordinates τ and σ, with the latter extending from σ = 0 (the endpoint of the string at r = ∞) to σ = σ max (the endpoint of the string at r = l 2 Λ QCD ). The variation of the action gives We have used the Wong equation (4.4) in order to obtain the above form. Also is the spacetime energy-momentum current carried by the string, D is the world-sheet covariant derivative built from the induced world-sheet metric h. Furthermore, n α is the normal to the world-sheet boundary; explicitly it reads (4.10) In particular, in the gauge (4.5), this becomes The string equation of motion which follows from the first line of (4.8) is The end of the string at σ = 0 is fixed at the Minkowski boundary r = ∞. Therefore, in a generic gauge, the coordinate r(τ, σ) of the string satisfies the Dirichlet boundary condition r(τ, σ = 0) = ∞ (4.13) for all times τ . Notice that this boundary condition is indeed compatible with the second line of (4.8) giving the variation of the world-sheet action. In the particular gauge (4.5), r is no longer a dynamical variable, as it becomes the world-sheet spatial coordinate. In that gauge, Eq. (4.13) is not a boundary condition anymore, but it merely specifies one of the string endpoints. Still at σ = 0, the other coordinates x µ satisfy the boundary condition −h τ τ n α P α µ + gQ a dx ν dτ F a µν = 0. (4.14) Once again, this is compatible with the second line of (4.8). Physically, Eq. (4.14) expresses the generalization of the equation of motion of a charged particle under the influence of the non-Abelian Lorentz force generated by the glasma fields. In our set-up, this force acts at the endpoint of the string and injects energy and momentum into the string in the Minkowski directions. Notice also that the string endpoint lives in a non-trivial space-time metric, cf. Eq. (2.21), due to the boundary sources representing hard fields 15 . We furthermore need to set boundary conditions at σ = σ max . For the radial coordinate r, we should clearly require r(σ = σ max , τ ) = l 2 Λ QCD , (4.15) which is indeed consistent with the variation of the action (cf. the third line in (4.8)). (In the gauge (4.5), this boundary condition on r should be thought of as a mere specification of the respective endpoint.) Regarding the other coordinates x µ , it is useful to recall the finding in previous studies [65][66][67][68][69], that the respective boundary conditions should not be set at the apparent horizon (which in our notations correspond to σ = σ max ), but rather at the world-sheet horizon. As known, a world-sheet horizon is bound to form whenever there is acceleration and/or an event horizon in the embedding space-time [70][71][72]. In our case, both situations may occur in general, as an event horizon dynamically emerges in AdS and moreover the endpoint of the string on the boundary can be accelerated by the 15 Note that G (s) µν blows up at the boundary r = ∞ like r 2 , however P r µ remains finite because ∂σx i (τ, r) (or ∂rx i (τ, r) in the gauge (4.5)) vanishes at the boundary with the right power of r. This means that the string always ends perpendicularly at the boundary. This is a familiar feature for trailing string solutions [65,66]. non-Abelian Lorentz force. In general the world-sheet horizon is located above the event horizon (i.e. closer to the boundary at r = ∞), but it may also coincide with the latter in special situations (e.g. for a string at rest). The emergence of a world-sheet horizon manifests itself as a singularity in the solution to the Nambu-Goto equation. Whenever this happens, one must choose the integration constant to ensure that this horizon is a regular point on the string, and not an endpoint (see [65,66] for an explicit realization of this idea in the case of the trailing string). Notice that the piece of the string which is located below the world-sheet horizon is causally disconnected from the string endpoint at the boundary. Accordingly, the energy and momentum fluxes across this world-sheet horizon should be interpreted as the energy and momentum losses by the heavy quark via radiation of soft and strongly-coupled quanta [69,[71][72][73][74][75].
To fully specify the solution to the Nambu-Goto equation (4.12), we also need initial conditions, say, at t = t 0 . For all the string points except for its endpoint at σ = 0, this initial condition refers to both the position and the velocity of the string points in the bulk. At σ = 0, however, it is sufficient to specify the position x i (t 0 ) of that endpoint: indeed, the boundary condition (4.14) determines the respective velocity including at the initial time, provided the initial condition for the color vector Q a is also specified.
In this respect, our set-up differs from other approaches, which employ the Dirichlet boundary condition at σ = 0, meaning that one needs to specify the entire trajectory x i (t) of the string endpoint on the boundary of AdS 5 , at any t. For instance, in Refs. [65,66], the string endpoint on the boundary is assumed to move at a constant velocity, under the combined action of its interactions with the plasma and of an external force. In that case, P r µ is calculated from the string profile, to obtain the (velocity-dependent) force acting on the string, which is then interpreted as (minus) the drag force acting on the heavy quark from the medium. In our approach, the endpoint of the string is fully dynamical and the energy loss (in particular the stopping distance) is obtained by solving the equations of motion.
In our particular physics set-up, where the soft sector is dynamically created at t > 0 by the hard sources, it is convenient to chose the initial conditions for the string at the collision time t = 0. A simple choice is to assume that, at t = 0, we have a vertical string extending from r = ∞ down to r = l 2 Λ QCD , and moving with a given velocity which is uniform along the string. This corresponds to an heavy quark with a specific initial kinetic energy in the center-of-mass frame of the collision.
To summarize, the string equations of motion (4.12) and the Wong equation Eq. (4.4) for the non-Abelian charge attached to the end-point at the boundary of AdS, supplemented with appropriate initial conditions, the boundary conditions (4.13), (4.14), and with the regularity condition at the world-sheet horizon as discussed above, should uniquely specify the dynamics of the string -particularly the energy-momentum loss towards the soft, strongly-coupled, sector. The string solution will also provide the trajectory x µ (t) of the heavy quark (meaning the trajectory of the end-point at the boundary of AdS) determining the time evolution of the non-Abelian charge current (4.7), which then can be used to calculate the additional energy loss towards the hard modes. In the future, we would like to develop an explicit algorithm for numerically solving the dynamics of this semi-holographic jet. The equivalent Polyakov action formulation as employed in [76] might be helpful in this regard.

Perspectives and future directions
In principle, our semi-holographic model can be used to calculate a wealth of experimental observables. In general, though, the associate numerical simulations may turn out to be rather cumbersome in practice. It would be therefore useful to dispose of a simplified picture, particularly for the late stages of the dynamics, where a state of local thermal equilibrium has been presumably reached and a hydrodynamical effective description becomes appropriate. It is indeed well known that, as far as the holographic part is concerned, the late time dynamics can be captured by an asymptotic hydrodynamic expansion with arbitrarily large number of transport coefficients [18,[77][78][79]. On the other hand the classical Yang-Mills equations are known to reduce to transport equations of the Boltzmann type as soon as the occupation numbers for the hard gluon modes become much smaller than 1/α s (which eventually occurs due to the longitudinal expansion) [80][81][82]. By appropriately combining these observations, one should be able to simplify our semi-holographic model in the final stages of the expansion, before the onset of hadronization.
The Boltzmann equation for the hard gluons emerging from our generalized Yang-Mills equations should be more general than the standard respective equation, in the sense of including additional contributions to the collision term which are associated with the exchange of bulk hydrodynamic modes via hard-soft couplings. Similarly, the hydrodynamics of the soft sector should be considered in a different effective metric that is modified as in (2.5), where the t µν of the hard-sector will be obtained from the Boltzmann equation. Thus we will obtain a novel and relatively simple description involving Boltzmann equations for the gluons coupled with the hydrodynamics of the fluid composed of the soft degrees of freedom, and living in an effective metric that is in turn determined by the gluon distribution which solves the Boltzmann equation. The only input of holography will be to provide the transport coefficients of the fluid composed of the soft degrees of freedom. Thus we reduce the problem of solving classical gravity equations to that of solving hydrodynamics living in a non-trivial background. At the same time, the classical Yang-Mills equations are replaced by a Boltzmann equation with a modified collision term. We hope that such a simplified picture will be sufficient to study observables like the collective flow.
Recently it has been shown that a generalization of the thermal fluctuation-dissipation relation [83,84] also holds for non-equilibrium holographic states which relax to static equilibrium and are sufficiently close to it -in the sense that one can describe the state by using either a hydrodynamic expansion in powers of derivatives, or an expansion in powers of the amplitude of the perturbation [85][86][87]. Some recent results indicate that such a generalized fluctuation-dissipation relation may also hold for the holographic system approaching a boost-invariant perfect fluid expansion and can be used to compute e.g. prompt photon and dilepton productions from the non-equilibrium QGP [88][89][90]. We hope to be able to design more direct tests of this state-independent non-equilibrium fluctuation-dissipation relation for the holographic sector at late time, by correlating suitable observables.
To conclude, more work is needed both on the theory side, to construct self-consistent numerical solutions to the coupled hard-soft equations and develop simplified versions of the model which are suitable at late time; and on the phenomenology side, to identify and compute new observables which are relevant for the physics of heavy-ion collisions. We intend to address such tasks in future work.