Leading order, next-to-leading order, and non-perturbative parton collision kernels: effects in static and evolving media

Energetic partons traveling in a strongly interacting medium lose energy by emitting radiation and through collisions with medium constituents. Non-perturbative, next-to-leading order and leading order collision kernels are implemented within AMY-McGill formalism. The resulting gluon emission rates are then evaluated and compared by considering scattering occurring in a brick of quark-gluon plasma, as well as in a realistic simulation of Pb-Pb, collisions at $\sqrt{s}=2.76$ ATeV using MARTINI. We find that the variations in quenching of hard partons resulting from using different kernels can be important, depending on the overall value of the strong coupling constant $\alpha_s$.


I. INTRODUCTION
The study of high-energy collisions of nuclei -"heavyions" -is the only practical way to study strongly interacting matter under extreme conditions of temperature and density, in terrestrial laboratories under controlled conditions. This active line of research seeks to understand quantum chromodynamics (QCD, the theory of the nuclear strong interaction) in all generality, and to answer fundamental questions about the phase diagram of QCD, and possibly to obtain precious information about the matter constituting the core of neutron stars, and the state of the universe a few microseconds after the Big Bang [1].
A vibrant global research program is under way, and the flagship facilities at the experimental energy frontier are RHIC (the Relativistic Heavy-Ion Collider, at Brookhaven National Laboratory) and the LHC (the Large Hadron Collider, at CERN). Several experiments performed at those laboratories have now revealed an exotic state of matter: the quark-gluon plasma (QGP) [2]. In the QGP, the partons constituting nuclear matter are free to roam over distances commonly considered macroscopic in hadronic physics. A considerable theoretical effort is currently being devoted to the detailed understanding of this QGP. In this context, one of the theoretical breakthroughs of the last couple of decades has been the realization that the collective behavior of the QGP can be modeled very faithfully by viscous relativistic fluid dynamics [3]. Hydrodynamics represents an effective theory of long wavelengths excitations, which can be extended to include nonequilibrium aspects with the introduction of transport coefficients such as shear (η) and bulk (ζ) viscosities. The connection between the phenomenology of relativistic heavy-ion collisions, fluid dynamics, and QCD is provided by the equation of state which is evaluated nonperturbatively in lattice QCD [4].
In addition to the measurements of one-body quantities like single-particle spectra of hadrons and of their collective motion -through the extraction of flow coefficients [3] -there exists a class of tomographic probes.
Those are penetrating and are able to reveal the inner workings of the QGP. Among those are QCD jets which form early in the heavy-ion collision, traverse and interact with the strongly interacting medium, fragment and are reconstructed by measuring their final states [5]. The very existence of QGP was revealed by measurements that highlighted the jet-medium interaction, through the analyses of di-hadron correlations [6,7]. The measurements of jets, of how they lose energy to the medium, and the theoretical modeling of this rich dynamics is therefore a central theme in modern high-energy nuclear physics .
The physics of the quark-gluon plasma is entering an era of precision science, where experiments and theory are precise and mature enough to provide quantitative characterization of this exotic material. One the lessons learned over the past decade or so is that a realistic dynamical modeling of the collisions system is key to the quantitative characterization of hot and dense strongly interacting matter. The current state-of-theart approaches combine partonic initial states, hydrodynamical evolution, and hadronic transport, while experimental final states are being connected with model parameters using techniques rooted in Bayesian statistics [34][35][36][37][38][39][40]. Those elaborate hybrid approaches model the background against which the tomographic signals evolve, and whose properties they probe. Precisely because QGP studies are maturing, this work focuses on one specific element of this chain: it aims to focus on the finite temperature jet-medium interaction as calculated within the Arnold-Moore-Yaffe (AMY) framework, and to explore quantitatively the empirical consequences of a recent extension of its partonic scattering kernel beyond the perturbative treatment [41,42].
The plan of this paper is the following: in Sec. II we discuss briefly the leading order radiative energy-loss rates and AMY-McGill formalism [28,43] as implemented in martini [44], as well as the new next-to-leading (NLO) order and nonperturbative (NP) scattering kernels introduced above in the context of a "brick test." In Sec. III we embed the new rates in a realistic music hydrody-namic simulation of Pb + Pb collisions at 2.76 ATeV. Finally in Sec. IV we discuss the results and their implications.

II. RADIATIVE ENERGY LOSS AND THE BRICK
A. Leading order rates The Arnold-Moore-Yaffe (AMY) formalism was proposed to calculate the gluon and photon emission from a thermalized plasma [45][46][47]. The aim was to compute the full leading order (in the strong coupling g s ) results for gluon (or photon) spectrum from a plasma, which meant the inclusion of bremsstrahlung and inelastic pairannihilation (or 1 → 2 channels). These channels are superficially at a higher order in α s than the Born contribution, but they receive parametric enhancements from the collinear singularities which are regulated by effective thermal masses. The end result of which is that these processes are found to be of the same order as the 2 → 2 channels. The inclusion of bremsstrahlung and inelastic pair-annihilation then requires the treatment of the Landau-Pomeranchuk-Migdal (LPM) effect, which is done with dynamic (i.e., nonstatic) scatterers. The complete leading order (LO) rates have been implemented in martini [44] along with the energy-loss associated with elastic channels [48]. martini then solves the Fokker-Planck-type rate equations for the evolving distribution of energetic partons, and contain a gain and a loss term: where p is the momentum of the original energetic parton and k is that of the emitted particle. In this form, the LO scattering rates are given by [44,49,50] dΓ (p, k) dk = g 2 s 16π p 7 where z ≡ k/p is the momentum fraction of the radiated particle, 1 1±e −k/T and 1 1±e −(p−k)/T are respectively the statistical factor for emitted parton and the remnant, which take into account the Bose enhancement (−) or Pauli blocking (+) effects, h ≡ (p × k) × p measures how collinear the emitted particle is to the mother, and is taken to be parametrically of O g s T 2 . The curly brackets contain the DGLAP parton splitting functions. Note that F (h, p, k) is the solution to the following integral equation where δE is the energy difference between the initial and final states and m i are the effective thermal masses induced by interactions with the medium. C s is the quadratic Casimir operator which for quarks(antiquarks) and gluons takes the values C s = C F = 4/3 and C s = C A = 3, respectively.
In the case where g → q +q, the prefactor (C s − C A /2) should multiply the term containing F (h − p q ⊥ ). Finally, C (q ⊥ ) is the transverse momentum broadening kernel: the differential rate of transverse (to the hard parton) momentum exchange with the medium. Assuming that the medium is constituted of dynamical charges, its value at leading order in α s is given by [51] where N c is the number of colors in the theory and N f is the number of fermion flavors under consideration. n B and n F are respectively the Bose and Fermi distribution, while p and p = p + q 2 ⊥ +2q ⊥ ·p 2(p−pz) are respectively the momenta of the medium particle before and after the scattering.

B. Jet-medium interactions at higher orders
More recently, efforts have been made to move beyond the LO calculation of rates and collisional kernels. The NLO scattering kernel has been computed by Caron-Huot [52] and used in a calculation of thermal photons [53] as well as for jet energy-loss in a medium within the AMY framework [54]. Strongly coupled media, however, may demand going beyond perturbation theory. In what concerns numerical calculations on the QCD lattice it had been known that going to a dimensionally reduced effective theory of QCD, known as EQCD or electrostatic QCD [55], enabled systematic higher-order 20 40 z 0 g → g + g (×100) 20 40 q → q + g (×260) calculations. Therefore, nonperturbative (NP) calculations of the elastic-scattering kernels via lattice simulations of EQCD have recently been performed, allowing for NP QCD input to jet-medium interactions and energloss. The NP scattering kernel has been computed on the lattice [56] for EQCD and matched to (3+1)D QCD for use in jet-medium energy-loss calculations [41] for an infinite medium. Medium-induced splitting rates were then computed using the new NP kernel and compared with to NLO and LO rates. Importantly, significant differences were observed between the NP and LO results [42].
In this work, we use the NLO and NP kernels developed in Ref. [41] in Eq. (3) and solve for the function F (h, p, k). The results of those computations are then implemented in martini together with our own recalculation of the LO rates [47], to quantify energy loss and to evaluate the consequences of the different scattering kernels. In the exploratory study performed here we use LO thermal masses, gluon emissions are still treated as collinear to the radiating parton, and collisional energy loss is treated at LO in the strong coupling [44].

C. Brick test of rates with LO, NLO and NP Collision Kernels
A first step in comparing the different scattering kernels is assessing the energy-loss rate. As is known from elementary statistics, very useful general information about a given distribution function may be obtained by studying the behavior of its moments. In this spirit, we consider the integrated rates and some of its moments where z ≡ k/p is the ratio between the momentum of the emitted or absorbed particle and that of the original energetic parton 1 . Then, z 0 g 4 s T ∆t and z 1 g 4 s T ∆t are respectively the difference between the number of parton emissions and that of parton absorptions, and the net energy-loss fraction divided by incident energy ("energy-loss fraction" from here on) within time interval ∆t. Also, z 2 z 0 − z 1 z 0 2 gives the variance of energy loss ratio per event, and Γ is the rate, computed using 1 The emission and absorption occur when z > 0 or z < 0, respectively.
dpT (GeV In what follows we refer to rates calculated using the NLO and NP kernels embedded in Eq. (2) as NLO or NP rates. Figure 1 compares the moments of the rates for q → q + g, g → g + g and g → q +q with the strong cou-pling g s and a factor of temperature scaled out for the LO, NLO and NP scattering kernels. At a given value of the scaled momentum, the magnitude of the net scattering rate evaluated with the NP kernel sits between that obtained with the LO and NLO expressions, clearly suggesting a slow convergence of the perturbative expansion. The slow convergence of many quantities evaluated in the framework of finite-temperature QCD is a common occurrence [1,52,53]. Depending on the value of the specific momentum (momentum divided by temperature) and on the moment, the results obtained with the LO, NLO, and NP kernels can differ by as much as 100%, as seen in Fig. 1 Observing the different moment distributions, the similarity of shapes obtained with the different kernels suggests a rescaling procedure. We therefore re-scale the curves by a constant factor, adjusted to match the ultraviolet behavior of the LO rate, and obtain the results shown as dashed curves in Fig. 1. It can be observed that the resulting scaled NP and NLO curves are very close to each other and match neatly with the high-momentum tail of the LO order rate. This is to be expected given the construction of NLO and NP kernels [54,56]. The major difference between the scaled rates and the LO re-sult is for p ≈ 10 T , where NLO and NP results can now be 10% (gluon to quark splitting) to 40% (gluon emission from gluon) higher than the LO rate, judging by the z 1 panels (i.e., the net energy-loss fraction). One can interpret this rescaling of the rates as using different values of α s , scaled appropriately relative to each other for the different rate sets. In fact, it is a common practice in phenomenological studies to tune α s according some observable that measures the medium modification of energetic partons such as charge hadron modification factor (R h ± AA ). Importantly, the fact that NLO and NP rates predict more quenching than the LO one at low transverse momentum is of phenomenological interest, and therefore of consequence. It is known [57] that simulations taking LO AMY rates with fixed α s predicted a p T -slope of R h ± AA that is more flat than the experimental results. Calculations with NLO and NP scattering kernels enhance the energy-loss for low momentum partons, which open the possibility of better agreement with experimental data.
To see the differences in a more detailed way, we perform a "brick test," using martini with the three rate sets. The brick is an idealized static QGP medium of infinite extent at a constant temperature. A quark propagates in the brick, and one studies how the medium influences the parton energy. Realistic nuclear collisions often create finite media with density and flow gradients: the brick dispenses with those complications in order to establish a baseline scenario. In what follows, unless stated otherwise, the traveling parton is a quark. We start with a parton of p = (50, 0, 0) GeV traveling within a QGP brick of temperature T = 0.3 GeV taking snapshots of its evolving spectrum. This representation of the time evolution enables an appreciation of how quickly the medium can transform the original δ function at a transverse momentum of 50 GeV. Figure 2(a) demonstrates the comparison of the rate sets when evaluated at the same value of strong coupling, α s = 0.3, a value appropriate for phenomenology at the BNL Relativistic Heavy Ion Collider (RHIC) and at the CERN Large HAdron Collider (LHC). The pattern of Fig. 1 where the NLO rate has a higher radiation rate relative to LO and NP is visible here, particularly for longer travel times where the quark jet subjected to NLO rates is clearly melting faster in the plasma. Figure 2 (b) implements the same conditions as Fig. 2 (a) but now allowing for elastic scattering. The inclusion of collisional energy-loss brings the distributions closer to each other without altering the observed hierarchy of energy-loss shown in Fig. 2 (a).
Scaling the value of α s for the brick test, as done before with the integrated rates in Fig. 1 according to makes the NP and NLO curves nearly indistinguishable from each other in Fig. 2 (c). After the scaling, the LO curve is well separated from the other two, particularly around the original jet momentum. Inclusion of the elastic-scattering channels for the case where the NLO and NP rates use a scaled value for α s , as can be seen in Fig. 2 (d), removes most of the difference between the LO rate and NLO, NP rates. There are some differences that remain for longer time spectra, in particular the NLO and NP curves slightly separating from each other while LO curve comes closer to both. While this was expected from Fig. 1 where the largest difference between the scaled (NLO, NP) rates and the LO result is expected to be for the smaller values of p/T , the rates plotted in that figure still imply some difference for p ≈ 100 T , even with a scaled α s . We can push the parton jet energy to the mini-jet domain while keeping the temperature constant and thus probe the region of the rates where the differences are more pronounced. Figure 3 is the equivalent calculation for a 10 GeV parton travelling through the same medium. Unlike the more energetic parton before, the difference between the rates is more clearly visible and they persist even after the inclusion of collisional energy-loss and rescaling of the coupling. This further reinforces the observation that the major difference between the rates lies in the infrared region.
We can further study the effect of the changing of the rates on the energy-loss of a energetic parton through a brick by considering the fractional energy-loss where p i is the initial momentum of the primary parton. Figure 4 provides a comparison of the un-scaled (α s = 0.3) and scaled NLO and NP curves to the LO result using the fractional energy-loss as defined in Eq. (8). Without collisional energy-loss we can clearly see that as the energy of the primary parton is reduced or temperature is increased the NLO and NP curves separate more from the LO result even when scaled. The closest the (NLO, NP) results get to the LO curve is for p = 100 GeV and T = 0.2 GeV (p/T = 500) which is in agreement with Fig. 1. Allowing for elastic scatterings as well as radiation in Fig. 5, while introducing a small difference between the NLO and NP curves, also reduces the difference between the new rates and the LO result. We can conclude that while for a fixed value of α s the NLO curve predicts more quenching relative to NP and LO rates, after scaling α s for each rate set we observe slightly more quenching for the NP kernel relative to others at lower values of p. This is marginally more pronounced when elastic scatterings are turned on.

III. ENERGY LOSS AND REALISTIC EVOLUTION
We now employ and study the new rates in a calculation of charged hadron nuclear modification factor, in an event-by-event simulation using realistic hydrodynamic backgrounds. N bin is the number of binary collisions in a nucleus-nucleus collision. We focus on 0 − 5% most central Pb + Pb collisions at beam energy 2.76 ATeV and compare with the corresponding experimental data. The hydrodynamic background used in this study is provided by music [58], a (3 + 1)D relativistic viscous hydrodynamic simulator of heavy-ion collisions. music concurrently solves the conservation equation of energymomentum tensor and the relaxation equations of shear stress tensor and bulk viscous pressure [59,60]. In this work, the fluid dynamics initial condition, at τ = 0.4 fm, of the stress-tensor is provided by two-dimensional ip-glasma calculation [61,62], and then the hydrodynamic quantities are evolved taking the HotQCD Collaboration equation of state [4], a constant shear viscosity to entropy density η/s = 0.13, and the same parametrization for the bulk viscosity as in Ref. [63].
The medium is then calculated in an event-by-event basis using music and stored as input for martini. Hard partons are evolved by martini through the medium and can lose (or gain) energy via their interactions with it. In this work we do not account for the medium response to the passage of energetic particles. The initial hard parton distribution is determined by pythia 8.209 [64] using the CMS underlying event tune (tune # 15) [65]. For our baseline pp calculations we use the CTEQ 6L1 [66] proton PDF while the Pb + Pb results are computed using nCTEQ15-np [67] nuclear PDFs. All PDF sets are provided by the LHAPDF 6.1.6 [68] package. While the hard scatterings are generated by pythia, their location is sampled in the transverse plane using the initial energy density as provided by ip-glasma. pythia then performs a full vacuum shower and passes the particle list to martini for evolution in the medium. It should be noted that the partons free-stream until the QGP medium is "formed" around them at τ = τ 0 which for us is taken to be τ 0 = 0.4 fm/c, the start time of the hydro. Evolution of momentum for a given parton is stopped if it can not be distinguished from a medium particle -defined here by the condition p < 4 T -or it is not in the fireball. Finally, energy-loss of jets in the hadron-cascade stage is not taken into account in this work. While preliminary studies of energy-loss in this stage indicate potentially significant contribution for hadron transverse momenta between 2 − 10 GeV [69] (outside of the p T range considered here), the interaction rates are much lower due to lower temperatures and therefore subdominant to energy-loss in the QGP stage.
For the sake of keeping this study simple, there are features we do not include in this work and leave for future consideration: correlated background subtraction (recoil effects) and coherence (formation time). Recoil effects concern the inclusion of medium particles that receive energy from energetic partons in an elastic-scattering event and get "promoted" to a jet particle. These, then, leave a hole in the medium which has to be subtracted at analysis level from the hadron and jet spectra. The main contribution of including recoil effects is in the better description of the jet transverse profile, an observable which we do not consider here. Their effect on the charged hadron R AA is restricted to the low p T region of the spectrum, below what we study here [57]. Formation time effect is time taken by the medium to distinguish the radiated parton from the incoming jet. This effect shows itself in the time dependence of radiative rates [70]. Finally, the different α s appearing in this study do not run. After all partons are frozen out, they are hadronized using the Lund string model as provided by pythia.  We present a sample of our baseline proton-proton calculation in Fig. 6. As can be seen, we have good agreement with the experimental measurements in three pseudorapidity windows, with data falling within 25% of the computed curves (lower panel). We then tune α s separately for all three rate sets, in order to fit the experimental data of R AA by minimizing  the chi-squared function, where the sum over i is over the experimental data points for p T > 10 GeV and the sum over s adds up all uncertainties. For the fitting procedure we use the experimental R AA in most central 0 − 5% 2.76 TeV Pb + Pb collisions from ALICE [73], ATLAS [71] and CMS [72].
With the fit results listed in Table. I, we note that the ratios between the fitted values of the strong-coupling are quite similar to the ratios predicted in Eq. (7), It confirms the understanding that difference between the three rate sets can be mostly resolved by rescaling the coupling constants. Meanwhile, the fact that they are slightly above the predicted values reveals the effect of elastic scattering, which is treated in the same way for all three rate sets. Fig. 7 (b) presents the results using fitted values for α s in each rate set. With a fitted α s the NLO and NP curves give nearly identical results to the LO results. The expected enhancement of p T -slope was not observed except the slight, although systematic, difference for p T < 30 GeV. There are two reasons for this. First, as shown in Fig. 1, one would expect a difference for partons with momentum p 30 T ∼ 10 GeV. The fragmented hadron from this parton is not in the region of interest in phenomenological studies. Second, the difference shown in Fig. 1 is the case only if there are no elastic scatterings present. In the presence of elastic scatterings, we could already see that the three rate sets would give very similar results (see Figs. 2(c),2(d) and 3(c) and 3(d)). It is known that partons with lower momentum receive significant modification from elastic channels and it is this factor that washes away the expected difference between the rate sets at low p T in charged hadron R AA .

IV. CONCLUSION AND OUTLOOK
The observation of energy-loss of jets and charged hadrons from heavy-ion collisions relative to a protonproton baseline remains an important signal of the creation of QGP. The dominant mode of energy-loss for energetic light-flavor quarks and gluons is gluon radiation as a result of soft scatterings with the medium. The interaction rate relies on the value of the strong-coupling as well as the collision kernel encoding information about momentum exchanges with the medium. It was found that the energy-loss by the QCD jets calculated using either the LO, NLO, or NP scattering kernels differed significantly when evaluated using the same value of α s . This is consistent with behavior of the radiative rates computed in Fig. 1. For example, the NLP curve for gluon emission from a gluon is ≈ 60% above the LO result while the NLO curve is nearly double the LO rate for p ∼ 10 T . This is an important conclusion of our study: for a given value of the strong coupling, the intrinsic differences between results computed with the different scattering kernels considered here can be large by any objective measure.
However, in phenomenological analyses of relativistic heavy-ion collisions, the overall value of α s is typically treated as a fitting parameter. Scaling the coupling value so that the new rates match the LO rate at large values of temperature-scaled momentum (p/T ) brings them much closer to the LO result, with NLO and NP nearly collapsing on each other and stand at only ≈ 40% above the LO curve for g → g + g at p ∼ 10 T . The brick test using the new radiative rates, excluding elastic scatterings with the medium, further demonstrated this effect with the evolving parton distribution from the NLO and NP radiative rates becoming nearly indistinguishable after the employed α s is rescaled (see Figs. 2(a),2(c) and 3(a), 3(c)). Including elastic channels [Figs. 2(d) and 3(d)] brings all three rate sets even closer to each other.
To study further the effect of the new collision kernels, we perform a realistic calculation of parton energy-loss using the new rates. This was done for the most centralt 0−5% Pb + Pb collisions at √ s = 2.76 ATeV using martini. The simulations, when run with equal and fixed values of α s , demonstrate the same behavior as the QGP brick test cases with the largest quenching done by the NLO rates. However, when we fit the value of the strong coupling, separately for three rate sets, to the charged hadron R AA results of ALICE, ATLAS and CMS in the same collision system, the resulting curves collapse on each other, becoming indistinguishable. This is because of the similarity of the scaled rates for large values of p/T as well as the effect of elastic scatterings in reducing the differences between different radiative rates. Moreover we observe basically no change in p T dependence of the charged hadron nuclear modification factor from using NLO and NP collision kernels in martini. An ad-ditional word of caution: martini as many jet-medium simulations, implements radiation as perfectly collinear with the emitting parton. This is an appropriate treatment for the LO kernel with a numerically small coupling (with opening angle of O (g s )), the transverse momentum dynamics will need to be scrutinized for all cases with g s 1.
As mentioned previously, we ignored the running of the strong-coupling in this study to avoid the associated complications and difficulties in interpreting our results. In a realistic run, however, α s will have a scale dependence with the physical consequence that the plasma becomes more "transparent" to a jet with high transverse momentum as compared with a simulation with a fixed α s [57]. Study of the new rate sets while accounting for a running coupling is currently underway.
Other means of distinguishing the behavior of the three kernels and their effect on the evolving distribution of hard partons in a plasma are analyses of multimessenger signals [74], or with more differential measurements such as those involving correlations. For the former, electromagnetic probes can offer potential access to jet-medium interactions [75], owing to the absence of final-state interactions with the medium. Work along those lines is in progress.