Multi-particle correlations and KNO scaling in the medium-induced jet evolution

We study the gluon distribution produced via successive medium-induced branchings by an energetic jet propagating through a weakly-coupled quark-gluon plasma. We show that under suitable approximations the evolution of the jet can be described as a classical stochastic process, which is exactly solvable. For this process, we construct exact analytic solutions for all the n-point correlation functions (the n-body densities in the space of energy). The corresponding results for the one-point and the two-point functions were already known, but those for the higher-point functions are new. These results demonstrate strong correlations associated with the existence of common ancestors in the branching process. By integrating these n-point functions over the gluon energies, we deduce the mean gluon multiplicity $\langle N\rangle$ as well as the higher moments $\langle N^p\rangle$ with $p\ge 2$. We find that the multiplicities of the soft gluons are parametrically large and show a remarkable regularity, known as Koba-Nielsen-Olesen (KNO) scaling: the reduced moments $\langle N^p\rangle/\langle N\rangle^p$ are pure numbers, independent of any of the physical parameters of the problem. We recognize a special negative binomial distribution which is characterized by large statistical fluctuations. These predictions can be tested in Pb+Pb collisions at the LHC, via event-by-event measurements of the di-jet asymmetry.

C Computation of h p (l) in two limiting cases 27 C.1 h p (l) in the limit π(pl) 2 1 27 C.2 h p (l) in the limit πl 2 1 28

Introduction
Motivated by extensive experimental studies of the energy loss by jets or leading hadrons in ultrarelativistic nucleus-nucleus collisions at RHIC and the LHC, and notably by the remarkable phenomenon known as 'di-jet asymmetry' [1][2][3][4][5][6][7][8][9], there has been a renewal of the interest in the theory and phenomenology of jet quenching -a concept which covers the ensemble of the modifications in the properties of a jet or of an energetic particle resulting from its interactions with a dense QCD medium, like a quark-gluon plasma. A substantial part of the recent developments refers to the evolution of a jet via multiple gluon emissions, as triggered by the collisions between the constituents of the jet and those of the medium. This is particularly interesting since, as observed in [10], the medium-induced gluon branchings have the potential to explain the striking pattern of the 'di-jet asymmetry', namely the fact that most of the energy lost by the subleading jet is taken away by many soft hadrons propagating at large angles w.r.t. to the jet axis.
Within perturbative QCD at weak coupling, the medium-induced jet evolution can be described as a classical stochastic process in which successive branchings are quasi-independent from each other [10][11][12][13] (see also Refs. [14][15][16][17][18] for earlier, related studies). In the most interesting physical regime, where the gluon radiation is triggered by multiple soft scattering, the branching rate is given by the BDMPSZ mechanism which takes into account the coherence between successive scatterings which contribute to a single emission (a.k.a. the Landau-Pomeranchuk-Migdal effect) [19][20][21][22][23][24][25][26][27][28][29]. Unlike the rate from bremsstrahlung in the vacuum, which depends only upon the splitting fraction z and favors soft (z 1) splittings, the BDMPSZ rate also depends upon the energy ω of the parent gluon and is such that it favors quasi-democratic branchings -that is, 1 → 2 branching processes for which the splitting fractions of the daughter gluons are comparable with each other: z ∼ 1 − z.
More precisely, for a jet propagating through the medium along a distance L, the BDMPSZ mechanism introduces a characteristic medium scale, the branching energy ω br (L) =ᾱ 2 sq L 2 , withq the jet quenching parameter (the rate for transverse momentum broadening via elastic collisions). The relatively soft gluons with energies ω ω br (L) have a probability of order one to disappear via democratic branchings and thus transmit their whole energy -via a mechanism akin to wave turbulence [10] -to a large number of very soft quanta, which are easily deviated to large angles via rescattering in the medium. In particular, if the medium is a quark-gluon plasma in thermal equilibrium, then the softest quanta produced by the branching process are expected to thermalize [30]. The energy taken away by these soft quanta can be identified with the energy lost by the jet towards the medium. In the experimental conditions at the LHC, the medium scale ω br (L) is expected to be relatively hard (of the order of a few GeV), yet significantly softer than the original energy E ≥ 100 GeV of the 'leading particle' (the parton initiating the jet). In this high-energy regime at E ω br (L), that will represent our main focus in that follows, the energy lost via soft quanta propagating at large angles is controlled by wave turbulence and is proportional to ω br (L) [10].
The average picture of the medium-induced jet evolution has been studied in great detail (at least, under suitable approximations to be later specified) [10,[30][31][32][33][34][35][36]. In particular, Ref. [10] presented exact, analytic, solutions for the gluon spectrum D(ω) = ω(dN/dω) and for the average energy lost by the jet at large angles. But this stochastic process is also expected to develop event-by-event fluctuations -say, in the number and the energy distribution of the branching products -, which have been less studied so far and which might be interesting too for the phenomenology.
The importance of fluctuations for a parton cascade generated via 1 → 2 gluon branchings is demonstrated by our experience with a jet evolving in the vacuum. In that case, the branchings are fully controlled (via the DGLAP dynamics) by the virtualities of the partons in the jet. The statistical properties of that process -that is, the mean gluon multiplicity and its higher momentshave been explicitly computed in the double-logarithmic approximation [37]. One has thus discovered [38] (see also Chapter 5 in the textbook [37]) the existence of large statistical fluctuations together with a remarkable regularity known as KNO scaling (from Koba, Nielsen, and Olesen [39]): in the large-virtuality limit, where the parton multiplicities are high, all the higher moments N p with p ≥ 2 are entirely fixed by the average multiplicity N . More precisely, the reduced moments N p / N p are pure numbers, independent of the virtuality, explicitly known in the approximations of interest.
Returning to our actual problem, that of a jet which propagates inside a dense medium, one may expect the associated fluctuations to be even larger, due to the 'democratic' nature of the mediuminduced gluon branchings: after each such a splitting, one loses trace of the parent parton, hence the ensuing cascade looks even more 'disordered' than for a jet which evolves in the vacuum. This will be confirmed and substantiated by our subsequent findings in this paper.
On the experimental side, the extensive studies of di-jet asymmetry at the LHC revealed that the 'missing energy' is balanced by a rather large number (10 to 15) of relatively soft hadrons (with p T between 0.5 and 2 GeV), which propagate in the hemisphere of the subleading jet, at large angles w.r.t. the jet axis [2,8,9]. Not surprisingly, this number shows large fluctuations event-by-event. However it seems difficult to also measure the multi-hadron correlations (like the p-body moments N p with p ≥ 2), due to the large background associated with the underlying event in a Pb+Pb collision.
On the theory side, there are only a couple of recent analyses of the importance of statistical fluctuations for the in-medium evolution of a jet (in particular, for the di-jet asymmetry) [40,41]. Ref. [40] presented a numerical study based on the Monte-Carlo event generator JEWEL [42], with very interesting conclusions: the di-jet asymmetry A J in central Pb+Pb collisions appears to be controlled by fluctuations in the branching process and not by the difference between the in-medium path lengths, L 1 and L 2 , of the two jets. In fact, the typical di-jet events generated by JEWEL are such that L 1 and L 2 are rather close to each other, a situation which according to the usual wisdom should lead to small values for A J . In spite of that, the numerical results exhibit a rather large fraction of events with large values for A J , including for the class of events where one enforces the condition L 1 = L 2 . This clearly demonstrates the importance of fluctuations.
In an independent, fully analytic, study which appeared soon after [41], we have for the first time computed the dispersion in the energy lost by the jet at large angles and in the multiplicity of soft gluons produced via medium-induced multiple branching. To that aim, we relied on an exact result for the gluon pair density, that we obtained under the same assumptions as used in previous studies of the average picture [10,[31][32][33][34][35][36]. Our results demonstrate that the fluctuations are huge: for both quantities alluded to above, the dispersion is parametrically as large as the respective mean value. In particular, if E denotes the energy lost by the jet event by event, then in the high-energy regime at E ω br , we found that σ E ∼ E ∼ ω br , where σ 2 E ≡ E 2 − E 2 . (We recall that ω br (L) =ᾱ 2 sq L 2 is the characteristic medium scale for multiple branching and E is the initial energy of the jet.) This in turn implies that the fluctuations in the medium-induced branching process can contribute to the di-jet asymmetry 1 at the same level as the difference in path-lengths between the two jets: (1.1) The first term in the r.h.s., proportional to the difference L 2 1 − L 2 2 , is the average contribution E 2 − E 1 and vanishes when L 1 = L 2 , as expected. But the second term, originating from fluctuations, is always nonzero and it dominates over the average piece whenever L 1 and L 2 are close to each other.
Clearly, the findings in Refs. [40,41] are consistent with each other, the reinforce each other and may together lead to a shift of paradigm concerning the physical origin of the di-jet asymmetry.
In this paper, we shall complete the study, started in Ref. [41], of multi-particle correlations in the gluon distribution produced via multiple medium-induced branchings. Our main new result is a set of exact, analytic, expressions for the p-body densities N (p) (x 1 , · · · , x p ) which describe the gluon distribution in the space of energy. Here, x i ≡ ω i /E is the energy fraction of one of the p measured gluons w.r.t. the leading particle. The respective expressions for p = 1 [the gluon density N (x), which is related to the spectrum D(x) via D(x) = xN (x)] and p = 2 [the gluon pair density N (2) (x 1 , x 2 )] were already known, as previously mentioned, but those for the higher-point correlations with p ≥ 3 are new. The 3-gluon density N (3) is shown in Eq. (3.7) and the generic p-body function N (p) in Eqs. (3.13)-(3.14). By inspection of these explicit results and of their derivation, one uncovers generic features and structural properties which shed more light on the physical picture for the in-medium jet evolution. For instance, all these correlations exhibit a 'geometric scaling' property: they depend upon the physical parameters α s ,q, L, and E only via the dimensionless ratio ω br (L)/E. The emerging physical picture for a typical event can be summarized as follows 2 : the jet is structured as an ensemble of 'mini-jets', i.e. gluon cascades generated via successive democratic branchings by 'primary partons' (gluons directly emitted by the leading particle) with relatively low energies ω ≤ ω br (L). Harder primary emissions, with ω br (L) ω < E, are possible as well, but they occur with a low probability (i.e. only in rare events) and do not give rise to mini-jets (since hard gluons cannot suffer democratic branchings). All the partons from a mini-jet are strongly correlated with each other, as they have a common ancestor to which they are linked via a succession of democratic branchings. On the other hand, different mini-jets are uncorrelated with each other, since the successive emissions of soft primary gluons are quasi-independent (indeed, the constraint of energy conservation plays only a minor role for the soft emissions).
By integrating the multi-gluon densities N (p) (x 1 , · · · , x p ) over their energy arguments x i , above some suitable infrared cutoff x 0 ≡ ω 0 /E, we shall deduce the gluon multiplicities -the average number N (ω 0 ) of gluons with energies ω ≥ ω 0 together with its higher moments 3 N p (ω 0 ). The lower cutoff ω 0 plays the role of an energy resolution scale. Without such a cutoff, the gluon multiplicity would be divergent, due to the rapid rise in the emission probability as ω → 0. On more physical grounds, one should observe that the 'ideal' (or 'turbulent') branching picture under consideration holds only for sufficiently high energies ω T , with T the average p T of the medium constituents (say, the temperature if the medium is a quark-gluon plasma). Hence, a physically meaningful cutoff satisfies ω 0 ∼ T ω br (L) E. With this choice, the multiplicities are parametrically large and dominated by the softest measured gluons, those with energies ω ∼ ω 0 . This last feature allows for relatively simple analytic approximations (see Sect. 4.1 for details).
Specifically, we shall find that N ∼ (ω br /ω 0 ) 1/2 1 and N p ∼ N p . This in particular implies KNO scaling: the reduced moments κ (p) ≡ N p / N p are pure numbers, independent of all the physical parameters of the problem, i.e. α s ,q, L, and E (see Eq. (4.6)). This feature is also shown by a jet evolving in the vacuum [37], but the respective distributions are significantly different: the statistical fluctuations are considerably larger for the medium-induced evolution. A precise way to characterize this difference is by comparing the corresponding probability distributions P(N ), which in the KNO regime are fully specified by the set of numbers κ (p) with p ≥ 1. For the medium-induced evolution, we shall recognize P(N ) as a specific negative binomial distribution (NBD) [43], that with parameter r = 2 (see Sect. 4.3 for details). This distribution is indeed broader (in the sense of developing larger fluctuations) than the one generated by a jet in the vacuum; the latter is known only numerically and can be viewed as an interpolation between the two NBD's with r = 2 and r = 3, respectively [37]. 2 We consider here the high-energy regime at E ω br (L), which is the most interesting one for the phenomenology at the LHC. For the corresponding picture at lower energies, E ω br (L), see the discussion in Sect. 3.
3 More precisely, we shall compute the factorial moments N (N −1) . . . (N −p+1) , but this distinction is unimportant in the high-multiplicity regime of interest, where N p.
The emergence of a NBD in relation with the medium-induced jet evolution is perhaps surprising, at the same level as the existence of exact, analytic, solutions for all the multi-gluon correlations. Notice indeed that we have a better analytic control for the evolution of the jet in the medium than in the vacuum, albeit the latter is a priori supposed to be a simpler problem. Without having a fully convincing explanation in that sense, we believe that both features could be related to the physics of wave turbulence, more precisely, to the existence of fixed-point solutions to the evolution equations for the multi-gluon densities N (p) (x 1 , · · · , x p ). On one hand, these fixed points greatly facilitate the search for analytic solutions. On the other hand, they determine the multi-gluon spectra at low energies x 1 and hence, in particular, the multiplicities of soft gluons. This paper is organized as follows. In Sect. 2 we briefly review the theoretical description of the medium-induced gluon branching as a Markovian stochastic process. In this context, we shall present for the first time the evolution equation obeyed by the p-body density N (p) (x 1 , · · · , x p ) with generic p ≥ 1. More details on the construction of this equation are presented in Appendix A. In Sect. 3 we present and discuss our exact results for the multi-gluon correlations N (p) (x 1 , · · · , x p ). We first recall the known results for p = 1 and p = 2 (but our physical discussion of the pair density N (2)

Master equations
We consider the parton cascade -to be subsequently referred as 'the jet', for brevity -which is generated via multiple gluon branchings by an incoming parton -the 'leading particle' (LP) -with initial energy E which crosses the medium along a distance (or time) L (the 'medium size'). We assume the LP to be on-shell at the time when it enters the medium, so that all the subsequent branchings are induced by its interactions, and the interactions of its descendants, with the constituents of the medium. We furthermore assume the medium to be a weakly coupled quark-gluon plasma in thermal equilibrium with temperature T E. The most important gluon emissions for what follows are those with intermediate energies, within the range T ω qL 2 , for which the formation times are much smaller than the medium size L, but much larger than the average mean free path between two successive collisions in the plasma. The rate for such medium-induced emissions can be computed in the multiple soft scattering approximation, with a result known as the BDMPSZ spectrum for a single gluon emission [19][20][21][22][23]. Moreover, as shown in [11][12][13], successive emissions can be treated as independent from each other, because the typical duration between two emissions (the 'branching time' to be introduced in Eq. (2.1) below) is parametrically larger than the formation time for individual emissions and, moreover, the coherence between the daughter gluons is efficiently washed out by the scattering in the medium. As a result, the jet evolution via medium-induced gluon branching can be described as a Markovian stochastic process [12,41].
Our goal throughout this paper is to study the energy distribution generated by this stochastic process, including fluctuations and correlations. To characterize this distribution, we shall compute the factorial moment densities N (p) (x 1 , x 2 , . . . , x p |τ ) for any p ≥ 1 and for generic values for the energy fractions x i ≡ ω i /E and for the 'reduced time' τ . The factorial moment density N (p) is roughly speaking the p-body density in the space of energy; this will be more precisely defined in Eq. (2.4) below. The reduced time is defined as τ ≡ t/t br (E), where t ≤ L is the actual time (or distance) travelled by the leading particle across the medium and the reference scale t br (E) is the 'democratic branching time' for the LP -that is, the 'lifetime' of the LP until it disappears via a quasi-democratic branching. By 'quasi-democratic' we mean a 1 → 2 gluon branching where the daughter gluons carry comparable fractions of the energy of their parent parton. For a parent gluon with energy ω, one has (see e.g. the discussion in [41]) withᾱ s ≡ α s N c /π (α s is the QCD coupling, assumed to be fixed, and N c is the number of colors) and q the 'jet quenching parameter' (the transport coefficient for transverse momentum diffusion).
We have anticipated here that, to the approximations of interest, the gluon distribution produced by the medium-induced jet evolution shows an interesting, geometric scaling, property: for given values of the energy fractions x i , the factorial moment densities N (p) (x 1 , x 2 , . . . , x p |τ ) depend upon the physical parameters of the problem -the travelled distance L, the transport coefficientq, and the original energy E of the LP -via a single, dimensionless, variable: the reduced time τ =ᾱ s L q/E. Accordingly, the gluon correlations are not modified when simultaneously changing, say, the medium size L and the energy E, but in such a way to keep constant the ratio L 2 /E.
The jet evolution via medium-induced multiple branching is a stochastic process whose dynamics is most economically expressed in terms of the generating functional with P n (x 1 , x 2 , · · · , x n |τ ) the probability density for having a state with n gluons with energy fractions x i (i = 1, . . . , n), at (reduced) time τ and u(x) an arbitrary 'source' function with support at 0 ≤ x ≤ 1. At τ = 0, we have just the LP, hence P n (τ = 0) = δ n1 δ(x 1 − 1). Probability conservation requires Z τ [u = 1] = 1 for any τ ≥ 0. The expectation value of an arbitrary observable is computed as where O n ≡ O(x 1 , x 2 , · · · , x n ) denotes the value of O in a particular event with n gluons. Strictly speaking, this probabilistic description requires an infrared cutoff (say, a lower limit on x), playing the role of an energy resolution scale, below which gluons cannot be resolved anymore. Indeed, the branching dynamics produces an infinite number of arbitrarily soft gluons and the 'state with exactly n gluons' is not well defined without such a cutoff. Any explicit construction of such a state, say via Monte-Carlo simulations, must involve an infrared cutoff on x, to be viewed as a part of the 'state' definition. On the other hand, the correlation functions of interest, like the p-body densities N (p) (x 1 , x 2 , . . . , x p |τ ), are insensitive to the unobserved, soft, gluons, hence they are independent of this cutoff. So long as one is solely interested in such semi-inclusive quantities, one can formally proceed without introducing any infrared cutoff.
Given the generating functional (2.2), the p-body densities of interest are obtained via functional differentiation w.r.t. the source field u(x) : In particular, for p = 1 one finds the gluon density in x-space, or gluon spectrum, and for p = 2, the density of pairs of gluons (with each pair being counted twice): The function N (p) (x 1 , x 2 , . . . , x p |τ ) with p ≥ 2 is totally symmetric under the permutations of the variables x i . The factorial cumulant densities, obtained from the logarithm of the generating functional, will be useful too, as they measure genuine correlations in the gluon distribution. In particular, The time evolution of all the p-body correlations is succinctly described by a single, functional, evolution equation for Z τ [u], which reads [41] where the kernel (the 'reduced' version of the BDMPSZ spectrum) is the differential probability per unit (reduced) time and per unit z for the splitting of a parent gluon with energy fraction x into a pair of daughter gluons with energy fractions zx and (1 − z)x, with 0 < z < 1 (the splitting fraction). The functional derivative δZ τ [u]/δu(x) in the r.h.s. of Eq. (2.8) plays the role of an 'annihilation operator' (it reduces by one the number of factors of u). The term quadratic in u within the square brackets describes the 'gain' effect associated with the branching process x → (zx, (1 − z)x), whereas the negative term linear in u is the corresponding 'loss' effect. By taking p functional derivatives in Eq. (2.8) and evaluating the result at u(x) = 1, it is straightforward to obtain the evolution equation obeyed by the p-th factorial moment. The respective equations for p = 1 and p = 2 have been presented in previous publications [10,41], but will be repeated here, for more clarity (and also because our present conventions are slightly different). They read and respectively Eq. (2.10) is homogeneous and must be solved with the initial condition N (x, τ = 0) = δ(x − 1). The first term in its r.h.s. describes the gain in the number of gluons at x due to emissions from gluons with x = x/z > x, whereas the second term describes the loss via the decay into softer gluons.
Eq. (2.11) is inhomogeneous; its r.h.s. involves the source term which describes the simultaneous creation of a pair of gluons with energy fractions x 1 and x 2 via the branching of a parent gluon with energy fraction x 1 + x 2 (with x 1 + x 2 ≤ 1 of course). After this splitting, the two daughter gluons evolve independently from each other and create their own gluon distributions (this evolution is described by the first line of Eq. (2.11)). Accordingly, the solution N (x, τ ) to Eq. (2.10) acts as a Green's function for Eq. (2.11): the solution to the latter with the initial condition N (2) (x 1 , x 2 |τ = 0) = 0 can be written as with a transparent physical interpretation: at some intermediate time τ , a gluon with energy fraction ξ 1 + ξ 2 splits into two gluons with energy fractions ξ 1 and respectively ξ 2 , whose subsequent evolutions generate two mini-jets which include the final gluons, x 1 and respectively x 2 . Note that the parent gluon with energy ξ 1 + ξ 2 is the last common ancestor (LCA) of the two measured gluons x 1 and x 2 . The equation obeyed by N (p) for generic p ≥ 1 will be derived in Appendix A and reads ∂ ∂τ 14) The similarity with Eq. (2.11) is quite manifest: the distribution N (p−1) for p − 1 particles acts as a source for the p-body density N (p) . Specifically, the source term in the above equation describes the simultaneous creation of the pair of particles x i and x j (taken among the p particles measured by N (p) ) via the splitting of one gluon with energy fraction x i + x j that was included in N (p−1) . Eq. (2.14) is formally solved by (with the initial condition N (p) (τ = 0)) where S (p) (ξ 1 , · · · , ξ p |τ ) denotes the source term in the r.h.s. of Eq. (2.14). Eq. (2.15) is truly a recursion formula, which expresses the p-body density N (p) in terms of the (p − 1)-th one.

Analytic solutions for the multi-gluon correlations
From now on, we shall focus on a slightly simplified version of the master equations introduced in the previous section, which has the virtue to allow for explicit, analytic, solutions, while at the same time keeping all salient features of the dynamics. This version is obtained by replacing the original kernel K(z) (i.e. the branching rate) with K 0 (z) ≡ 1/[z(1 − z)] 3/2 . The simplified kernel preserves the pole structure of the exact kernel at z = 0 and z = 1, hence it generates a very similar evolution. This is confirmed by numerical solutions to Eq. (2.10) using both forms of the kernel [33,35]. For this simplified kernel, one was able to obtain exact, analytic solutions for the gluon spectrum N (x, τ ) [10] and the pair density N (2) (x 1 , x 2 |τ ) [41]. In what follows, we shall first briefly review these known solutions and thus introduce a physical picture for the medium-induced jet evolution which will later be refined by our new results for the higher p-point functions.

The gluon spectrum
The gluon spectrum N (x, τ ) corresponding to the simplified kernel K 0 (z) reads [10] N For relatively small times πτ 2 1, i.e. t t br (E), this spectrum exhibits a pronounced peak near x = 1 which describes the leading particle together with a power tail N (x, τ ) τ /x 3/2 at x 1, which describes soft radiation. The shift 1 − x p ∼ πτ 2 in the position of the peak measures the typical energy lost by the LP via radiation, where the width of this peak δx p ∼ πτ 2 describes the 'broadening', i.e. the uncertainty in the energy of the LP due to fluctuations in the radiation process. In physical units, with the notations ∆E p ≡ E(1 − x p ) and δE p ≡ Eδx p , one finds, parametrically, This 'branching energy' ω br (t), which is independent of E and much smaller than it (remember that we are in the regime where πτ 2 1), is the characteristic energy scale for the onset of multiple branching : there is a probability of order one to emit a gluon with energy ∼ ω br (t) during a time interval t. In particular, a quasi-democratic branching occurs with probability of O(1) during t provided the energy ω br (t) of one (any) of the daughter gluons is comparable to the energy ω of their parent; this condition ω ∼ ω br (t) implies t ∼ t br (ω), with t br (ω) introduced in Eq. (2.1).
The estimates in Eq. (3.2) can be physically understood as follows [41]: during the relatively small time t t br (E), the LP cannot undergo a democratic branching, rather it radiates a large number of very soft gluons with energies ω ω br (t) together with a number of O(1) of harder gluons with ω ∼ ω br (t). The latter control the energy lost by the LP in a typical event, hence ∆E p (t) ∼ ω br (t). Besides, the fluctuations in the number of hard gluons are of O(1) as well (since successive hard emissions are quasi-independent), hence also the broadening δE p (t) must be of order ω br (t).
After being emitted by the LP, the primary gluons with energies ω ω br (t) are bound to undergo democratic branchings, because the corresponding branching times obey t br (ω) t < L. Via successive democratic branchings, they generate parton cascades ('mini-jets') and thus gradually transfer their energy to softer and softer quanta, and eventually to the medium. This mechanism for energy transfer is extremely efficient, since characterized by wave turbulence [10]: the rate for energy flow from one parton generation to the next one is independent of the generation (i.e. of x). This is visible in the special form of the gluon spectrum (3.1) at small x, namely N (x) ∝ x −3/2 : this power law represents a turbulent fixed point for the branching process. That is, the gain and loss terms in the r.h.s. of Eq. (2.10) mutually cancel for this particular spectrum, meaning that there is no net accumulation of energy at any intermediate values of x: after a time of order t br (ω), the whole energy ω that was initially carried by a primary gluon ends up into arbitrarily soft quanta (x → 0).
More precisely, this 'turbulent' branching picture applies so long as the gluon energies remain larger than the temperature of the medium: the very soft gluons with energies of order T can efficiently exchange energy and momentum with the medium constituents, via elastic collisions, and thus thermalize [30]. As a result, the energy ω of the primary gluon is ultimately transmitted to the medium, in the form of many soft quanta which emerge at large angles.
This discussion shows that the energy transfer from the LP to the medium via democratic cascades is quasi-deterministic. Accordingly, the statistics of the energy lost by the LP, as encoded in Eq. (3.2), also applies to the energy lost by the jet as a whole, via soft radiation at large angles: both the average energy loss and its dispersion are of order ω br (t), since controlled by the hardest 'primary' emissions which can occur in a typical event [41]. The above discussion applies so long as t t br (E), which is the relevant situation for the phenomenology of di-jet asymmetry at the LHC [1,2]: the initial energy E of the LP is so high (E ≥ 100 GeV) that the corresponding branching time t br (E) is much larger than the available medium size L 6 fm.
For completeness, let us also consider the situation at larger times t t br (E), or τ 1. This applies to jets with a smaller overall energy E and also to the 'mini-jets' generated by primary gluons, as previously discussed. When πτ 2 ∼ 1, the exponent in Eq. (3.1) becomes of O(1) for any value of x, hence the peak disappears from the spectrum: as expected, the LP undergoes a democratic branching. When further increasing τ , the support of the spectrum shrinks towards small values x 1/πτ 2 1, while its strength is rapidly suppressed, as e −πτ 2 . This confirms that all the gluons with x > 0 (more precisely, x > T /E) disappear via democratic branchings.

The gluon pair density
The gluon pair density has been computed in [41], with the following result (a streamlined derivation of this result will be given below, in Appendix B.1) To better appreciate the physical interpretation of this result, it is useful to have a look at its derivation. The difference of two exponentials in the r.h.s. has been generated via the following integral, where the integration variable τ has the same meaning as in Eq. (2.13): this is the splitting time for the last common ancestor (LCA). The integrand in Eq. (3.4) happens to be a total derivative w.r.t. τ , hence the result of the integration comes from the two endpoints. The first term in the r.h.s. is generated by the upper limit τ = τ and describes a process in the which the splitting of the LCA occurs very late, near the time of measurement; the respective exponential is recognized as the probability for the LCA (with energy fraction x 1 + x 2 ) to survive over a time τ (compare to Eq. (3.1)). The second exponential similarly refers to an early splitting (τ = 0) and can be interpreted as the probability that both daughter gluons survive until they are finally measured at time τ . For relatively large measurement times, πτ 2 1, both exponentials are small, yet the first process (late splitting) dominates over the first one, since it is less likely for two particles to survive than for a single one.
In practice though, we are more interested in the small-time regime τ 1, as appropriate for the phenomenology of jets at the LHC. When πτ 2 1, the pair density (3.3) develops a peak near x 1 + x 2 = 1, corresponding to the case where one of the measured gluons is the LP. (See the plot in Fig. 1 for an illustration.) But the most interesting situation is when both x 1 and x 2 are small, x 1 , x 2 1, as generally the case for radiation. In that case, Eq. (3.3) reduces to where we have also used the corresponding estimate for the gluon spectrum, that is, N (x, τ ) τ /x 3/2 (cf. the discussion after Eq. (3.1)). In spite of its factorized structure, the pair density in Eq. (3.5) does still encode strong correlations, as shown by the following argument: the genuine 2-particle correlation is measured by the cumulant pair density where the second estimate, valid in the regime of Eq. (3.5), is parametrically as large as N (2) (x 1 , x 2 |τ ) in that regime. As argued in [41], this correlation comes from processes where the LCA is itself soft (ξ 1 + ξ 2 1 in Eq. (2.13)). The other possibility, namely that the LCA be the LP, does not generate correlations, since successive, soft emissions by the LP are nearly independent from each other.
To leading order in perturbation theory -by which we mean the expansion in the number of gluon emissions (or, equivalently, the perturbative solution to the master equations obtained via iterations) -, the result in Eq. (3.5) receives contributions from the three processes shown in Fig. 2, each of them involving two soft emissions. It is quite easy to explicitly compute the respective contributions (this requires two iterations of Eqs. (2.10) and (2.11)) and thus check that they sum up to the result shown in Eq. (3.5); each channel contributes 1/3 of the total result and the net correlation comes from the leftmost channel, where the LCA is itself soft. However, in reality, there are arbitrarily many other processes, involving the emissions of unresolved gluons with energies ω ω br (t) (or energy fractions x τ 2 ), which contribute to the same accuracy. Indeed, as discussed after Eq. (3.2), the probability for such a soft emission is of order one, so one can include arbitrarily many of them without modifying the perturbative accuracy of a calculation. The fact that such additional emissions do not modify the net result beyond the 'naive' leading-order calculation is again a consequence of wave turbulencethat is, of the precise cancellation between 'gain' and 'loss' contributions to the r.h.s. of Eq. (2.11).

The 3-gluon correlation
In the Appendix B, we shall construct an inductive argument allowing one to compute the factorial moments N (p) for arbitrary p. Before discussing the general case, in the next subsection, let us here present and discuss the respective result for p = 3. This is conveniently written as where we have introduced the complimentary error function, The expansion of the error function for small values of its argument will be useful too for what follows: The special linear combination of error functions appearing in Eq. (3.7) has been generated via the following double time integral (compare to Eq. (3.4)) τ 0 dτ 2 (3.10) In this integral τ 1 represents the splitting time for the branching generating the 2-point correlation N (2) , which subsequently acts as a source for the 3-point function (cf. Eq. (2.14) with p = 3). Furthermore, τ 2 refers to the branching which creates the final 3-point correlation; that is, this is the same as the integration variable τ in Eq. (2.15). As emphasized in Eq. (3.10), the integrand can be written as a total derivative w.r.t. τ 1 ; hence the integral over τ 1 is trivial and the subsequent integral over τ 2 is recognized as the definition of the error function. The 3 terms within the square brackets in Eq. (3.7) correspond to the 3 possible combinations of late and early emissions. The first term, which yields the dominant contribution at large times πτ 2 1, represents processes where both splittings occur very late, close to the time of measurement: τ 1 τ 2 τ . Accordingly, this term is proportional to the survival probability for a common ancestor with energy fraction x 1 +x 2 +x 3 . [This becomes obvious after using the asymptotic behavior of erfc(a) at large a 1, cf. Eq. (3.8).] The last error function in Eq. (3.7) corresponds to the case where both emissions occur very early, τ 1 τ 2 0, while the intermediate one, to the case where τ 1 0 and τ 2 τ (an early emission plus a late one). These last 2 terms are strongly suppressed at late times, since proportional with the survival probabilities for systems of 3 and, respectively, 2 particles.
As already explained, the most interesting situation for the phenomenology at the LHC is the small-τ regime at πτ 2 1. In this regime and for small energy fractions 4 x i 1 (corresponding to the bulk of the radiation), one can use the expansion of the error function in Eq. (3.9) to find Note that the would-be dominant terms, linear in τ , have cancelled out between the various error functions, hence the dominant contribution is cubic in τ , as it should in order to be consistent with factorization. In spite of this factorized structure, the small-time result in Eq. (3.11) does still encodes strong correlations, as obvious when computing the respective cumulant (cf. Eq. (2.7)): (3.12) To lowest order in perturbation theory, the result in Eq. (3.11) receives contributions from processes involving the emission of exactly 3 soft gluons, as illustrated in Fig. 3. In particular, the net correlation in Eq. (3.12) is generated by the leftmost process, where the 3 measured gluons have a soft common ancestor. As already discussed in relation with the 2-point function, this result (3.11) is truly nonperturbative, in that it receives contributions from processes with arbitrarily many unresolved soft gluons. However, all contributions beyond order τ 3 exactly cancel because of the fine cancellations between gain and loss terms -i.e. because of wave turbulence.

The generic p-body density
The previous discussions of the gluon pair density N (2) , in Sect. 3.2, and of the 3-body density N (3) , in Sect. 3.3, were quite similar to each other and this similarity has inspired us an induction argument 4 The behavior near the kinematical limit at x1 + x2 + x3 = 1, where the argument τ / √ 1 − x1 − x2 − x3 of the error functions can be large even for small values of τ , is not that interesting since the 3-point function N (3) is strongly suppressed in that limit, as we shall see in Sect. 3.4. Figure 3. Processes contributing to the production of 3 soft gluons to lowest order in perturbation theory. The process on the left, where the LCA is itself soft, is the only one to generate a 3-particle correlation.
which allows for the calculation of the higher-point correlations. This argument will be presented in detail in the Appendix B. Here, we shall merely present the final result and discuss some physical consequences.
Specifically, our result for N (p) can be conveniently written as with the function h p ( ) defined by the following multiple integral, The time-dependence of N (p) is fully encoded in the function h p and hence it enters only via the scaling variable ≡ τ / 1 − p j=1 x j . Accordingly, the structure of the factorial moment N (p) is remarkably simple: this is essentially the product of the power-like spectrum produced by wave turbulence (meaning one factor of 1/x 3/2 i for each external leg), which controls the gluon distribution at small x (x i 1), times a scaling function which describes the time-dependence of the multi-gluon correlation and also its behavior near the kinematical limit at p j=1 x j = 1. We did not attempt to analytically perform the time integrations in Eq. (3.14) for generic values of p (the first such integral, over 1 , is of course trivial since the integrand in Eq. (3.14) is a total derivative). But, clearly, this multiple integral is well suited for numerical calculations and also for analytic approximations, as we shall now discuss.
We first consider the situation where the scaling variable is large, 1. This includes the large-time regime, τ 1, but also the behavior near the kinematical limit at p j=1 x j = 1 for generic values of τ . As already seen on the examples of the 2-point and 3-point functions and it is also intuitive by inspection of Eq. (3.14), the dominant behavior in this limit comes from processes where all the relevant splittings occur as late as possible: i for any i = 1, 2, . . . , p − 1. Indeed, such configurations minimize the exponent of the Gaussian within the integrand of Eq. (3.14). We thus expect an asymptotic behavior h p ( ) ∝ e −π 2 , proportional to the survival probability of the last common ancestor (with energy fraction x 1 + x 2 + · · · + x p ) over a time of order τ . This is confirmed by the manipulations in the Appendix C.2, which more precisely yield This holds up to corrections suppressed by inverse powers of and/or exponentials factors like e −4π 2 . This approximation yields This is in agreement with our previous results for the p-point functions with p = 1, 2, 3. Notice that, unlike the gluon spectrum (3.1) and the pair density (3.3), which exhibit a leading-particle peak in the vicinity of the kinematical limit, the higher-point correlations with p ≥ 3 do not show such a peak -rather, they rapidly vanish when approaching that limit. We now turn to the more interesting regime at small times τ 1 and soft produced particles x i 1. Clearly, this means τ 1. It is then tempting to evaluate the integrations in Eq. (3.14) by using the small-argument expansion of the exponential. For this to be justified for generic values of p (including larger values p 1), one however needs the stronger condition π(p ) 2 1. Under this stronger assumption, the dominant behavior is indeed obtained by replacing the Gaussian by unity. One thus finds (see the Appendix C.1 for details) 17) and therefore for πτ 2 1 p 2 and x i 1 , (3.18) in agreement with our respective results for p = 1, 2, 3. This power-like multi-particle spectrum, with the characteristic exponent 3/2, is the consequence of wave turbulence for the jet problem at hand. In spite of the factorized structure of Eq. (3.18), this result still encodes genuine correlations, as already discussed on the examples of the 2-point and 3-point functions. These correlations are generated by processes where all the p measured particles have a soft common ancestor -a primary gluon whose evolution via democratic branchings has generated the mini-jet to which all the measured particles belong. As explained, the genuine correlations can be isolated by computing the factorial cumulant.
As an example, we here show the corresponding result 5 for p = 4 : To illustrate the previous approximations for the function h p ( ), we display in Fig. 4 the function h 3 ( ) (whose explicit form can be easily inferred by comparing Eqs. (3.7) and (3.13)) together with its approximate versions at 1, cf. Eq. (3.15), and at 1, cf. Eq. (3.17). 5 One has (with compact notations, whose meaning should be obvious): C

Gluon multiplicities and KNO scaling
In the previous section, we have studied the detailed energy distribution of the medium-induced radiation, as characterized by the factorial moment densities N (p) (x 1 , · · · , x p |τ ). In what follows, we shall 'integrate out' the distribution in energy in order to deduce the statistics of the gluon multiplicities -the total number of gluons and its fluctuations.

Multiplicities for soft gluons
If one attempts to compute the average number of gluons by integrating the gluon spectrum (3.1) over x, that is, N (τ ) = 1 0 dx N (x, τ ), then one faces a severe infrared (x → 0) divergence, due to the strong, power-like, enhancement in the gluon density at small x : N (x, τ ) ∝ 1/x 3/2 . This argument shows that the total gluon number is not a meaningful observable, since the radiation produces infinitely many soft gluons. Rather, it makes sense to compute the total number of gluons with energies larger than some minimal value ('infrared cutoff') ω 0 , meaning with energy fractions x ≥ x 0 , where x 0 ≡ ω 0 /E. (We recall that E is the initial energy of the leading particle.) A natural value for ω 0 exists on physical grounds: this is the characteristic energy of the medium, say, its temperature T if the medium is a weakly-coupled quark-gluon plasma in (or near) thermal equilibrium.
Indeed, the ideal branching dynamics considered so far strictly applies only so long as the energies of the gluons from the cascade remain much larger than this medium scale, ω T . On the other hand, when ω ∼ T , the dynamics is modified, first, by the elastic collisions between the gluons from the jet and the constituents of the medium and, second, by the non-linear effects associated with the relatively large gluon occupation numbers (which for ω ∼ T become of O(1), since one cannot distinguish anymore between gluons from the jets and those from the medium). As a result of such modifications, the soft gluons with ω T are expected to thermalize, which in turn will stop the branching process (due to the detailed balance between splitting and recombination processes) [30].
For the kinematical conditions at the LHC, one has T ω br (L) E, where we recall that ω br (L) =ᾱ 2 sq L 2 is the characteristic energy scale for the onset of multiple branching (cf. the discussion after Eq. (3.2)). The first inequality ensures that one has a sufficiently large phase-space at soft momenta for the jet evolution via multiple branchings to be fully developed. The second inequality implies that the medium is relatively thin, L t br (E), so the leading particle survives in the final state -it radiates soft gluons with ω ω br (L), but it cannot undergo a democratic branching.
In practice, it is useful to chose the infrared cutoff ω 0 to be much smaller than ω br (L) (in order to probe the physics of multiple branchings) but still larger than T (to be able to distinguish the particles from the jet from those in the surrounding medium). Accordingly, we shall compute the p-order factorial moment of the multiplicity as follows, with the lower cutoff x 0 satisfying T /E < x 0 τ 2 1. The result is strongly sensitive to the precise value of x 0 , yet in what follows we shall identify observables which are independent of this cutoff. For the physical interpretation of the subsequent results, it is useful to recall the relations When x 0 τ 2 , i.e. ω 0 ω br (L), the multiple integral in Eq. (4.1) is controlled by its lower limit, i.e. x i ∼ x 0 for any i = 1, . . . , p, due to the copious production of soft gluons via multiple branching.
One can therefore estimate this integral by using the approximate version of N (p) (x 1 , · · · , x p |τ ) valid when x 1 + x 2 + · · · + x p 1, that is (cf. Eq. (3.13)), Then the integrals in Eq. (4.1) become trivial and yield As it should be clear from the above, this approximation properly counts the soft gluons produced via radiation, but it ignores the relatively hard gluons with τ 2 x ≤ 1 and notably the leading particle. Hence, this should be a good approximation so long as the gluon multiplicities are very high, in particular N (τ, x 0 ) 1, which is strictly true so long as x 0 τ 2 . In practice though, this remains a good approximation up to larger values x 0 ∼ τ 2 provided one pays attention not to include the LP when counting the multiplicity -a condition which is easy to fulfill in the experimental situation, where the LP can indeed be distinguished from its soft products of radiation so long as τ < 1. To illustrate this, we compare in Fig. 5 the prediction of Eq. (4.4) for the average multiplicity N (τ, x 0 ) against the exact respective result [41] from which we subtract the contribution of the LP. That is, Eq. (4.4) with p = 1 is compared to N (τ, x 0 ) − e −πτ 2 , with N (τ, x 0 ) given by Eq. (4.4) in Ref. [41]. As visible in this figure, the approximation in Eq. (4.4) is indeed accurate up to x 0 ∼ τ 2 .

KNO scaling
A remarkable feature of the result in Eq. (4.4) is the fact that it scales as a power of 1/x 0 , with an exponent proportional to p. This reflects the power-like spectrum of the factorial moments at small x i , cf. Eq. (4.3), which we recall is a (turbulent) fixed point of the evolution via multiple branching. In turn, this implies that the dependence upon x 0 cancels out when constructing the reduced moments, The reduced moments (4.5) exhibit geometric scaling : they depend upon the physical parameters L,q, and E only via the dimensionless variable τ =ᾱ s L q/E. This scaling takes a particularly simple form at sufficiently small τ , when one can use the approximation (3.17) for h p (τ ). In that case, the ratio h p (τ )/[h 1 (τ )] p (p + 1)/[2(p − 1)!] becomes independent of τ , hence the reduced moment (4.5) is a pure number, which depends only upon p : This property is known as KNO scaling (from Koba, Nielsen, and Olesen [39]). In a previous publication [41], we have obtained this result (4.6) for the particular case p = 2 and conjectured the emergence of KNO scaling for generic values of p. The present analysis confirms the existence of this scaling, clarifies the limits of its validity (cf. the inequalities in the r.h.s. of (4.6)), and also specifies the corresponding value for κ (p) for any p ≥ 2. As we shall shortly explain, this value for κ (p) is quite special and allows us to identify the probability distribution for gluon multiplicities at small τ . In practice, the KNO scaling is restricted to rather small values of τ -the more so, the larger the value of p. This limitation is indeed visible in Fig. 6, where we plot the ratio κ (p) (τ )/κ (p) (0) computed according to Eq. (4.5) as a function of τ , for p = 2 and p = 3. The respective KNO prediction, namely κ (p) (τ )/κ (p) (0) = 1 (cf. Eq. (4.6)), is seen to be satisfied only at very small values of τ . Besides, the deviation from it with increasing τ starts earlier (and grows faster) for p = 3 than for p = 2. Physically, this can be understood as follows: the scaling occurs so long as all the configurations (in the sense of branching trees within the parton cascade) that contribute to the simultaneous production of a set of p particles survive with a probability of order one. This requires e −π(pτ ) 2 ∼ O(1), hence τ 2 1/(πp 2 ). The KNO scaling should be better and better satisfied with increasing the jet energy E (since τ 2 ∝ 1/E), but it is unclear whether this can be observed within the current experimental conditions at the LHC.
Incidentally, under the assumptions of Eq. (4.6) it is possible to obtain a relatively compact expression for the factorial moments themselves (and not only for their ratios), namely As indicated by the above notations, the soft gluon multiplicities in this high-energy regime become independent upon the energy E of the LP, but only depend upon the (large) ratio ω br (L)/ω 0 between the characteristic medium scale ω br (L) and the energy resolution scale ω 0 . Moreover, the multiplicities are parametrically large in this regime, hence the factorial moments can be identified with the ordinary moments: one has N (N − 1) . . . (N − p + 1) N p for any fixed value of p. Let us finally emphasize a rather general and also robust conclusion of the previous discussion, namely, the fact that the fluctuations in the soft multiplicities are very large. This becomes more transparent if one considers the respective cumulants, like N (N − 1) − N 2 , which are a direct measure of correlations. Recall indeed that all the cumulants vanish for the Poisson distribution, which describes independent emissions. But for the present distribution, the cumulants are non-zero and large -in fact, as large as possible: indeed, they are parametrically as large as the factorial moments themselves. For instance, using Eq. (4.5) for x 0 τ 2 , one finds e πτ 2 2πτ 2 for τ 1 , where we have also used the approximations (3.17) and (3.15) for small and respectively large values of τ . This behavior is illustrated in Fig. 7. The result 1/2 at small τ is an immediate consequence of Eq. (3.6) and comes from processes where the 2 measured particles belong to a same mini-jet, i.e. they have a soft common ancestor. For larger τ 1, the second cumulant is even larger than the disconnected piece N 2 ; this is a consequence of the fact that the large-τ behavior of N (N − 1) is controlled by special configurations where the splitting of the last common ancestor occurs towards the end of the evolution (cf. the discussion after Eq. (3.4)). Similar conclusions apply for the higher cumulants with p > 2. They show that the correlations in the gluon distribution at small x are indeed very strong. A more precise characterization of these correlations at small τ (i.e. in the high energy limit) will be presented in the next subsection.

A negative binomial distribution
As previously mentioned, the soft gluon multiplicities in the small-τ (or high-energy) limit do not only exhibit KNO scaling, but they also suggest a remarkable probability distribution, that can be directly read off Eq. (4.6) : this is the negative binomial (or Pascal) distribution with parameter r = 2. To be more specific, let us first recall that the negative binomial distribution (NBD) involves 2 free parameters, the average particle numbern ≡ n and a positive integer 6 r whose meaning will be shortly explained. The associated probability law reads: One possible interpretation for the random variable n which is consistent with this distribution is as follows: n is the number of failures before the occurrence of a prescribed number r of successes in a sequence of independent Bernoulli trials (see e.g. [43,44] for details). Using (4.9) one can easily compute the associated generating functional, from which it is straightforward to deduce the factorial moments: For r = 2, this is in agreement with Eq. (4.6), as anticipated.
We are not aware of any fundamental physical reason for the emergence of this particular NBD in the jet problem at hand. Moreover, it should be clear from the above that the actual multiplicity distribution inside the jet is generally different from (and more complicated than) a NBD: in general, the reduced moments (4.5) depend upon an additional parameter τ , which is real and positive. They reduce to the simple form in Eq. (4.6) only for sufficiently small values of τ , where what we mean by 'sufficiently small' depends upon p -the precise condition becomes more and more restrictive with increasing p. This means that, strictly speaking, there is no fixed value of τ , even if arbitrarily small, for which the actual multiplicity distribution is truly equivalent to the NBD with r = 2.
Moreover, even in the physical regime where our approximations (4.6) and (4.7) make sense, they do not probe the details of the NBD for generic values of n andn, but only its tail at large multiplicities, n 1 andn 1. Indeed, Eq. (4.7) applies only in the high-multiplicity (and high-energy) regime at ω 0 ω br (L) E/p 2 . Vice-versa, the expression (4.6) for the reduced moment κ (p) is also generated by the simplified version of Eq. (4.9) (with r = 2) valid at large n and largen, that is, nP n (n, 2) 4ρ e −2ρ , ρ ≡ n n . (4.12) Notice that, in this approximation, the quantitynP n scales as a function of ρ = n/n, a property which is sometimes used as the definition of KNO scaling. Let us rapidly check that (4.12) implies indeed the result (4.6) for κ (p) in the limit where n 1. In this limit, one can ignore the difference between factorial and ordinary moments and replace the sum over n by an integral: (4.13) 6 Generalizations to real values of r are also possible, but they are not useful for our present purposes.
But even if limited, this relation between the distribution of soft gluons within the in-medium jet and the NBD is conceptually interesting, as we now explain. This becomes clearer when the present problem is compared to the evolution of a jet in the vacuum, as driven by its virtuality. In that case too, one found that the soft particle multiplicities obey KNO scaling, with reduced moments κ (p) that have been explicitly computed (within a double logarithmic approximation) [37]. In particular, one found κ (2) = 4/3, which via Eq. (4.11) appears to be consistent with a NBD with parameter r = 3. This identification is not fully right -for p ≥ 3, the respective values κ (p) start deviating 7 from those of the NBD with r = 3 -, yet this is representative for the statistical properties of the in-vacuum jet evolution. The smaller the value of r, the broader is a negative binomial distribution 8 . Hence, also by this argument, we conclude that the multiplicity distribution created by the medium-induced evolution is considerably broader -in the sense of developing larger statistical fluctuations -than that associated with a jet which propagates in the vacuum.

Conclusions
In this paper we have investigated the multi-particle correlations in the gluon distribution generated via medium-induced multiple branchings by an energetic jet propagating through a weakly coupled quark-gluon plasma. We have demonstrated that, under suitable approximations, the jet evolution can be described as a stochastic branching process which is exactly solvable: for this process, we have obtained exact analytic results for all the p-body gluon densities in the space of energy. The corresponding results for p = 1 (the gluon spectrum) [10] and p = 2 (the gluon pair density) [41] were already known in the literature, but those for the higher-point correlations with p ≥ 3 are new. By integrating these densities over the energies of the gluons, above an infrared cutoff ω 0 which plays the role of the resolution scale, we have deduced the factorial moments N (N − 1) . . . (N − p + 1) (L, ω 0 ) which characterize the distribution of the gluon multiplicity.
The results that we have thus obtained have interesting physical consequences, which could be observed in the experiments. They demonstrate large multiplicities for the soft gluons together with strong correlations associated with the existence of common ancestors. While such correlations were to be expected in the context of a branching process, they appear to be significantly stronger than for the corresponding process in the vacuum (the DGLAP evolution of a jet driven by the virtuality of the leading particle) [37]. This reflects the fundamental difference between the respective branching laws: unlike the rate for bremsstrahlung in the vacuum, which is scale invariant and favors soft splittings (i.e. splittings where the daughter gluon carries only a small fraction of the energy of its parent parton), the BDMPSZ rate for medium-induced gluon branching involves the dimensionful transport coefficient q and favors democratic splittings (at least for the sufficiently soft gluons -those whose energies are softer than the characteristic medium scale ω br =ᾱ 2 sq L 2 ). As a result, the medium-induced branchings are strongly biased towards soft energies. This leads to the abundant production of soft gluons with ω ω br . Moreover, any such a gluon becomes the seed of a 'mini-jet' produced via a sequence of democratic branchings. All the gluons within the same mini-jet are correlated with each other, as they have a common ancestor.
The expressions for the soft gluon multiplicities and the associated correlations become particularly suggestive in the high energy limit E ω br : they are independent of the energy E of the leading particle and scale as powers of the large ratio ω br /ω 0 ; specifically, N p ∝ ω br /ω 0 p/2 . This in particular implies that N p grows with the medium size as L p . It furthermore implies that the reduced moments κ (p) ≡ N p / N p are pure numbers, independent of any of the physical parameters of the problem. This property is known as KNO scaling. The specific value of κ (p) , cf. Eq. (4.6), is instructive too: it implies that the associated probability distribution is a special negative binomial distribution, which is over-dispersed -i.e., it features large fluctuations. This appears to be more dispersed than the distribution produced by the DGLAP evolution of a jet propagating in the vacuum. It would be interesting to search for confirmations of these results in the experimental data at the LHC -notably, in the distribution of soft particles at large angles in the context of di-jet asymmetry. Most likely, it should be difficult to see indications of the KNO scaling: on one hand, our prediction in that sense relies on idealized theoretical assumptions; on the other hand, it is notoriously difficult to experimentally measure multi-particle correlations in nucleus-nucleus collisions, due to the large background associated with the underlying event. Yet, some of the qualitative consequences of our results are already consistent with the LHC data for di-jet asymmetry: the fact that the multiplicities of soft hadrons propagating at large angles are large and characterized by large eventby-event fluctuations [1,2,8].
Let us finally recall the main assumptions underlying our present analysis, to clarify its limitations and suggest directions of improvement for further studies (which however will most likely require numerical simulations). First, the medium has been described as a static quark-gluon plasma, characterized by a homogeneous transport coefficientq. In view of the phenomenology, one must extend this set-up to an expanding medium, with a time-dependent and possibly also space-dependent distribution forq, that could be dynamically generated via elastic collisions. A suitable framework in that sense is provided by the AMY kinetic equations [45] which lie at the basis of the Monte-Carlo event generator MARTINI [46]. Second, we have limited ourselves to a leading-order formalism within perturbative QCD at weak coupling. However, as recently understood, there are important quantum corrections -notably, double-logarithmic corrections toq from medium-induced radiation [47][48][49][50][51] and thermal corrections of O(g) to the kinetic equations [52,53] -that are by now available and could be used to improve our current estimates. In particular, the double-logarithmic corrections are non-local and introduce an additional dependence upon the medium size L, in the form of an 'anomalous dimension' [48][49][50].
Furthermore, we have neglected the virtuality of the leading particle and the associated vacuumlike radiation. Whereas one does not expect the vacuum-like radiation to directly contribute to the energy loss at large angles, it may indirectly do so, by producing additional sources at small angles. Besides, this is itself a random process, which introduces additional fluctuations. We are not aware of any analytic formalism allowing for the simultaneous treatment of the parton virtualities and the in-medium collisions, and of their combined effect in triggering radiation. However, this becomes possible (at least, modulo some approximations) within Monte-Carlo event generators like JEWEL [42]. It was indeed in that Monte-Carlo context that the importance of fluctuations for the in-medium jet evolution has been first pointed out [40]. That numerical approach has also the virtue to allow for additional sources of fluctuations which are known to be important for the phenomenology, like those in the geometry of the hard process and the distance L travelled by a jet through the medium (a quantity that has been treated as fixed in our analysis). We hope that our present analytic findings, although obtained in a somewhat idealized set-up, will inspire more systematic Monte-Carlo studies and thus open the way to realistic applications to the phenomenology. Based on that, it is easy to guess that where A p some proportionality constant that only depends on p. We can check explicitly that if this assumption is fulfilled for p − 1 then it is also fulfilled for p with A p = p 2 A p−1 . Using this together with A 2 = 1/2π, we find where f has the property that f (l|l 1 , · · · , l i , · · · , l j , · · · , l p−1 ) = f (l|l 1 , · · · , l j , · · · , l i , · · · , l p−1 ) for any i and j. Therefore we can rewrite the expression for h p as We can also apply the change of variables l i → lλ i to obtain (C.3) Using this equation, we shall now study the two limiting cases of physical interest.
C.1 h p (l) in the limit π(pl) 2 1 The equation (C.3) can be simplified if one can substitute the exponential by 1 for all the values of λ i inside the domain. This will precisely happen when π(pl) 2 1 as this is the maximum value of the argument of the exponential, in this case we can use that (C.7) When πl 2 1, the integral will be dominated by the region y i ∼ 1 (C.9)