Langevin dynamics of heavy quarks in a soft-hard factorized approach

By utilizing a soft-hard factorized model, which combines a thermal perturbative description of soft scatterings and a perturbative QCD-based calculation for hard collisions, we study the energy and temperature dependence of the heavy quark diffusion coefficients in Langevin dynamics. The adjustable parameters are fixed from the comprehensive model-data comparison. We find that a small value of the spatial diffusion coefficient at transition temperature is preferred by data $2\pi TD_{s}(T_{c}) \simeq 6$. With the parameter-optimized model, we are able to describe simultaneously the prompt $D^{0}$ $R_{\rm AA}$ and $v_{2}$ data at $p_{\rm T}\le8$ GeV in Pb--Pb collisions at $\sqrt{s_{\rm NN}}=2.76$ and $\sqrt{s_{\rm NN}}=5.02$ TeV. We also make predictions for non-prompt $D^{0}$ meson for future experimental tests down to the low momentum region.

By utilizing a soft-hard factorized model, which combines a thermal perturbative description of soft scatterings and a perturbative QCD-based calculation for hard collisions, we study the energy and temperature dependence of the heavy quark diffusion coefficients in Langevin dynamics. The adjustable parameters are fixed from the comprehensive model-data comparison. We find that a small value of the spatial diffusion coefficient at transition temperature is preferred by data 2πT Ds(Tc) 6. With the parameter-optimized model, we are able to describe simultaneously the prompt D 0 RAA and v2 data at pT ≤ 8 GeV in Pb-Pb collisions at √ sNN = 2.76 and √ sNN = 5.02 TeV. We also make predictions for non-prompt D 0 meson for future experimental tests down to the low momentum region.

I. INTRODUCTION
Ultrarelativistic heavy-ion collisions provide a unique opportunity to create and investigate the properties of strongly interacting matter in extreme conditions of temperature and energy density, where the normal matter turns into a new form of nuclear matter, consisting of deconfined quarks and gluons, namely quark-gluon plasma (QGP [1]). Such collisions allow us to study the properties of the produced hot and dense partonic medium, which are important for our understanding of the properties of the universe in the first few milliseconds and the composition of the inner core of neutron stars [2][3][4]. Over the past two decades, the measurements with heavy-ion collisions have been carried at the Relativistic Heavy Ion Collider (RHIC) at BNL and the Large Hadron Collider (LHC) at CERN [5][6][7], to search and explore the fundamental properties of QGP, notably its transport coefficients related to the medium interaction of hard probes.
Heavy quarks (HQs), including charm and bottom, provide a unique insight into the microscopic properties of QGP [8][9][10][11][12][13]. Due to large mass, they are mainly produced in initial hard scatterings and then traverse the QGP and experience elastic and inelastic scatterings with its thermalized constituents [14,15]. The transport properties of HQ inside QGP are encoded in the HQ transport coefficients, which are expected to affect the distributions of the corresponding open heavy-flavor hadrons. The resulting experimental observables like the nuclear modification factor R AA and elliptic anisotropy v 2 of various D and B-mesons are therefore sensitive to the HQ transport coefficients, in particular their energy and temperature dependence. There now exist an extensive set * lish@ctgu.edu.cn † sunfei@ctgu.edu.cn ‡ xiewei@ctgu.edu.cn § xiongw@ctgu.edu.cn of such measurements, which allow a data-based extraction of these coefficients. In this work, we make such an attempt by using a soft-hard factorized model (see Sec. III) to calculate the diffusion and drag coefficients relevant for heavy quark Langevin dynamics.
A particularly important feature of the QGP transport coefficients is their momentum and temperature dependence, especially how they change within the temperature region accessed by the RHIC and LHC experiments. For instance the normalized jet transport coefficientq/T 3 was predicted to present a rapidly increasing behavior with decreasing temperature and develop a near-T c peak structure [16]. The subsequent studies [17][18][19][20][21][22] seem to confirm this scenario. Another important transport property, shear viscosity over entropy density ratio η/s, also presents a visible T -dependence with a considerable increase above T c [23]. Concerning the HQ diffusion and drag coefficients, there are also indications of nontrivial temperature dependence from both the phenomenological extractions [24][25][26][27][28][29][30] and theoretical calculations [20,21,[31][32][33][34]. The difference among the derived hybrid models is mainly induced by the treatment of the scale of QCD strong coupling, hadronization and non-perturbative effects [35,36]. See Ref. [37][38][39] for the recent comparisons.
The paper is organized as follows. In Sec. II we introduce the general setup of the employed Langevin dynamics. Section III is dedicated to the detailed calculation of the heavy quark transport coefficients with the factorization model. In Sec. IV we show systematic comparisons between modeling results and data and optimize the model parameters based on global χ 2 analysis. With the parameter-optimized model, the energy and temperature dependence of the heavy quark transport coefficients are presented in Sec. V, as well as the comparisons with data and other theoretical calculations. Section VI contains the summary and discussion.

II. LANGEVIN DYNAMICS
The classicial Langevin Transport Equation (LTE) of a single HQ reads [14] with i, k = 1, 2, 3. The first term on the right hand side of Eq. 1b represents the deterministic drag force, F i drag = −η D p i , which is given by the drag coefficient η D (E, T ) with the HQ energy E = p 2 + m 2 Q and the underlying medium temperature T . The second term denotes the stochastic thermal force, F i thermal = C ik ρ k / √ dt, which is described by the momentum argument of the covariance matrix C ik , together with a Gaussian-normal distributed random variable ρ k , resulting in the uncorrelated random momentum kicks between two different time scales κ T and κ L are the transverse and longitudinal momentum diffusion coefficients, respectively, which describe the momentum fluctuations in the direction that perpendicular (i.e. transverse) and parallel (i.e. longitudinal) to the propagation. Considering the Einstein relationship, which enforces the drag coefficient starting from the momentum diffusion coefficients as [40] The parameter ξ denotes the discretization scheme of the stochastic integral, which typically takes the values ξ = 0, 0.5, 1, representing the pre-point Ito, the mid-point Stratonovic, and the post-point discretization schemes, respectively; d = 3 indicates the spatial dimension. In the framework of LTE the HQ-medium interactions are conveniently encoded into three transport coefficients, i.e. η D , κ T and κ L (Eq. 4). All the problems are therefore reduced to the evaluation of κ T /L (E, T ), which will be mainly discussed in this work. Finally, we introduce a few detailed setups of the numerical implementation. The space-time evolution of the temperature field and the velocity field are needed to solve LTE (Eq. 1a and 1b). Following our previous analysis [27,41], they are obtained in a 3+1 dimensional viscous hydrodynamic calculation [42], with the local thermalization started at τ 0 = 0.6 fm/c, the shear viscosity η/s = 1/(4π) and the critical temperature T c = 165 MeV (see details in Ref. [42]). When the medium temperature drops below T c , heavy quark will hadronize into the heavy-flavor hadrons via a fragmentation-coalescence approach. The Braaten-like fragmentation functions are employed for both charm and bottom quarks [43,44]. An instantaneous approach is utilized to characterize the coalescence process for the formation of heavy-flavor mesons from the heavy and (anti-)light quark pairs. The relevant coalescence probability is quantified by the overlap integral of the Winger functions for the meson and partons, which are defined through a harmonic oscillator and the Gaussian wave-function [45], respectively. See Ref. [40] for more details.

III. MOMENTUM DIFFUSION COEFFICIENTS IN A SOFT-HARD FACTORIZED APPROACH
When propagating throughout the QGP, the HQ scattering off the gluons and (anti-)partons of the thermal deconfined medium, can be characterized as the two-body elementary processes, with p 1 = (E 1 , p 1 ) and p 2 are the four-momentum of the injected HQ (Q) and the incident medium partons i = q, g, respectively, while p 3 and p 4 are for the ones after scattering. Note that the medium partons are massless (m 2 = m 4 ∼ 0) in particular comparing with the massive HQ (m 1 = m 3 = m Q in a few times GeV). The corresponding four-momentum transfer is (ω, q ) = (ω, q T , q L ). The Mandelstam invariants read with q ≡ | q | E 1 for small momentum exchange. The transverse and longitudinal momentum diffusion coefficients can be determined by weighting the differential interaction rate with the squared transverse and longitudinal momentum transfer, respectively. It yields and As the momentum transfer vanishes (|t| → 0), the gluon propagator in the t-channel of the elastic process causes an infrared divergence in the squared amplitude |M 2 | t−channel ∝ 1/t 2 , which is usually regulated by a Debye screening mass, i.e. t → t − λm 2 D with an adjustable parameter λ [46,47]. Alternatively, it can be overcome by utilizing a soft-hard factorized approach [48,49], which starts with the assumptions that the medium is thermal and weakly coupled, and then the interactions between the heavy quarks and the medium can be computed in thermal perturbation theory. Finally, this approach allows to decompose the soft HQ-medium interactions with t > t * , from the hard ones with t < t * . For soft collisions the gluon propagator should be replaced by the hardthermal loop (HTL) propagator [50,51], while for hard collisions the hard gluon exchange is considered and the Born approximation is appropriate. Therefore, the final results of κ T /L include the contributions from both soft and hard components.
As discussed in the QED case [48], µ + γ → µ + γ, the complete calculation for the energy loss of the energic incident heavy-fermion is independent of the intermediate scale t * in high energy limit. However, in the QCD case, there is the complication that the challenge of the validity of the HTL scenario, m 2 D T 2 [49], due to the temperatures reached at RHIC and LHC energies. Consequently, in the QCD case, the soft-hard approach is in fact not independent of the intermediate scale t * [33,47]. In this analysis we have checked that the calculations for κ T /L are not sensitive to the choice of the artificial cutoff t * ∼ m 2 D . In the next parts of this section, we will focus on the energy and temperature dependence of the interaction rate Γ at leading order in g for the elastic process, as well as the momentum diffusion coefficients (Eq. 7 and 8) in soft and hard collisions, respectively.
A. κ T /L in soft region t * < t < 0 In soft collisions the exchanged four-momentum is soft, √ −t ∼ gT (λ mf p ∼ 1/g 2 T [52]), and the t-channel longwavelength gluons are screened by the mediums, thus, they feel the presence of the medium and require the resummation. Here we just show the final results, and the details are relegated to A. The transverse and longitudinal momentum diffusion coefficients can be expressed as respectively, with the HQ velocity v 1 = | p|/E 1 and the transverse and longitudinal parts of the HTL gluon spec-tral functions § are given by B. κ T /L in hard region tmin < t < t * In hard collisions the exchanged four-momentum is hard, , and the pQCD Born approximation is valid in this regime. In analogy with the previous part we give the κ T /L results directly, and the detailed aspects of the calculations can be found in B. The momentum diffusion coefficients reads The integrations limits and the short notations are shown in Eq. B8-B14.

C. Complete results in soft-hard scenario
Combining the soft and hard contributions to the momentum diffusion coefficients via The spectral function involves only the low frequency excitations, namely the Landau cut, while the quasiparticle excitations is irrelevant in this regime. See Eq.A8 for more details.
while κ sof t T /L (E ≡ E 1 , T ) is given by Eq. 9 and 10, and κ hard−Qi T /L is expressed in Eq. 13 and 14 for a given incident medium parton i = q, g. Adopting the post-point discretization scheme of the stochastic integral, i.e. ξ = 1 in Eq. 4, the drag coefficient η D (E, T ) can be obtained by inserting Eq. 15 into Eq. 4.

IV. DATA-BASED PARAMETER OPTIMIZATION
Following the strategies utilized in our previous work [40,53], the two key parameters in this study, the intermediate cutoff t * and the scale µ of running coupling (Eq. A14a), are tested within a wide range of possibility and drawn constrains by comparing the relevant charm meson data with model results. We calculate the corresponding final observable y for the desired species of D-meson. Then, a χ 2 analysis can be performed by comparing the model predictions with experimental data In the above σ i is the total uncertainty in data points, including the statistic and systematic components which are added in quadrature. n = N − 1 denotes the degree of freedom (d.o.f ) when there are N data points used in the comparison. In this study, we use an extensive set of LHC data in the range p T ≤ 8 GeV: D 0 , D + , D * + and D + s R AA data collected at mid-rapidity (|y| < 0.5) in the most central (0 − 10%) and semi-central (30 − 50%) Pb-Pb collisions at √ s NN = 2.76 TeV [54,55] and √ s NN = 5.02 TeV [56], as well as the v 2 data in semicentral (30 − 50%) collisions [57][58][59].
We scan a wide range of valus for (t * , µ): A total of 20 different combinations were computed and compared with the experimental data. The obtained results are summarized in Tab. I. The χ 2 values are computed separately for R AA and v 2 as well as for all data combined. To better visualize the results, we also show them in Fig. 1, with left panels for R AA analysis and right panels for v 2 analysis. In both panels, the y-axis labels the desired parameters, |t * | several models achieve χ 2 /d.o.f 1 with |t * | m 2 D 1.5 and widespread values of µ: 1 ≤ µ πT ≤ 3. This suggests that R AA appears to be more sensitive to the intermediate cutoff while insensitive to the scale of coupling constant. For the v 2 , it clearly shows a stronger sensitivity to µ, which seems to give a better description (χ 2 /d.o.f 1.5 − 2.0) of the data with µ πT = 1. It is interesting to see that R AA data is more powerful to constrain |t * | m 2 D , while v 2 data is more efficient to nail down µ πT . Taken all together, we can identify a particular model that outperforms others in describing both R AA and v 2 data simultaneously with χ 2 /d.o.f = 1.3. This one will be the parameteroptimized model in this work: |t * | = 1.5m 2 D and µ = πT .

V. RESULTS
In this section we will first examine the t * dependence of κ T /L (Eq. 15) for both charm and bottom quark. Then, the relevant energy and temperature dependence of κ T /L (E, T ) will be discussed with the optimized parameters. For the desired observables we will perform the comparisons with the results from lattice QCD at zero momentum limit, as well as the ones from experimental data in the low to intermediate p T region.  κ T and κ L are, as expected, identical at zero momentum limit (E = m Q ), while the latter one has a much stronger energy dependence at larger momentum. Furthermore, κ T /L (E, T ) behave a mild sensitivity to the intermediate cutoff t * [33]. Because the soft-hard approach is strictly speaking valid when the coupling is small m 2 D T 2 [49], thus, the above observations support the validity of this approach within the temperature regions even though the coupling is not small. In Fig. 3, charm quark κ T (left) and κ L (middle) are evaluated with the optimized parameters, |t * | = 1.5m 2 D and µ = πT , including both the soft (dotted blue curves) and hard contributions (dashed black curves), at fixed temperature T = 0.40 GeV (upper) and at fixed energy E = 10.0 GeV (lower). It is found that the soft components are significant at low energy/temperature, while they are compatible at larger values. The combined results (solid red curves) are presented as well for comparison. With the post-point scheme (ξ = 1 in Eq. 15), the drag coefficients (right) behave (1) a nonmonotonic dependence, in particular on the temperature, which is in part due to the improved treatment of the screening in soft collisions, and in part due to the procedure of inferring η D from κ T /L [15]; (2) a weak energy dependence at E 2m c = 3 GeV. Similar conclusions can be drawn for bottom quark, as shown in Fig. 4. Figure 5 presents the transport coefficient of charm quark,q = 2κ T , at fix momentum p = 10 GeV (solid red curve). It can be seen thatq/T 3 reaches the maximum near the critical temperature, and then followed by a decreasing trend with T , providing a good description of the light quark transport parameter (black circle points). The results from various phenomenological extractions and theoretical calculations, including a phenomenological fitting analysis with the Langevin-transport with Gluon Radiation (LGR; dotted blue curve [40]), a LO calculation with a Linearized Boltzmann Diffusion Model (LIDO ‡ ; dashed black curve [60]), a nonperturbative treatment with Quasi-Particle Model (QPM in Catania; dot-dashed green curve [61]), a novel confinement with semi-quark-gluon-monopole plasma approach (CUJET3; shadowed red band [18,20,21]), are displayed as well for comparison. Similar temperature dependence can be observed except the CUJET3 approach, which shows a strong enhancement near T c regime.
The scaled spatial diffusion coefficient describes the low energy interaction strength of HQ in medium [63], and it can be calculated by substituting Eq. 4 into Eq. 17. The obtained result for charm quark (m c = 1.5 GeV) is displayed as the solid red curve in Fig. 6. It is found that a relatively strong increase of 2πT D s (T ) from crossover temperature T c toward high temperature. Meanwhile, the v 2 data prefers a small value of 2πT D s near T c , 2πT D s (T c ) 3 − 6 § , which is close to the lattice QCD calculations [64][65][66][67][68]. The relevant results from other theoretical analyses, such as LGR [40], LIDO [60], Catania [61] and CUJET3 ¶ [18,20,21], show a similar trend but with much weaker temperature dependence.  lations are done with FONLL initial charm quark spectra and EPS09 NLO parametrization for the nPDF in Pb [27], and the pink band reflects the theoretical uncertainties coming from these inputs. It can be seen that, within the experimental uncertainties, the model calculations provide a very good description of the measured p T -dependent R AA data for various charm mesons. Concerning the results in Pb-Pb collisions at √ s NN = 5.02 TeV, as shown in Fig. 8, a good agreement is found between the model and the measurement at p T 6 GeV/c, while a slightly larger discrepancy observed at larger p T . ‡ LIDO results are shown with only elastic scattering channels. § This range is estimated from the testing models as displayed in the right panels of Fig. 1. ¶ CUJET3 results are obtained by performing the energy interpolation down to E = m Q .      (Color online) Transport coefficient,q/T 3 (T ), of charm quark from the various calculations, including: the soft-hard factorized approach (solid red curve), the LGR model with data optimized parameters (dotted blue curve [40]), a Bayesian ayalysis from LIDO (dashed black curve [60]), a quasi-particle model from Catania (dot-dashed green curve [62]), CUJET3 (shadowed red band [18,20,21]) and JET Collaboration (black circle points [19]) at p = 10 GeV.
B-hadron decays (dashed curve), reflecting a weaker inmedium energy loss effect of bottom quark (m b = 4.75 GeV), which has larger mass with respect to that of charm quark (m c = 1.5 GeV). We note that the future measurements performed for the non-prompt D 0 R AA and v 2 , are powerful in nailing down the varying ranges of the model parameters (see Sec. IV), which will largely improve and extend the current understanding of the inmedium effects. LGR (dotted blue curve [40], LIDO (dashed black curve [60]), Catania (dot-dashed green curve [61]), CUJET3 (shadowed red band [18,20,21]) and lattice QCD calculations (black circle [64], blue triangle [65], pink square [66], red inverted triangle [67] and green plus [68]). The result for D-meson (long dashed pink curve [69]) in the hadronic phase is shown for comparison.

VI. SUMMARY
In this work we have used a soft-hard factorized model to investigate the heavy quark momentum diffusion coefficients κ T /L in the quark-gluon plasma in a data-driven approach. In particular we've examined the validity of this scenario by systematically scanning a wide range of possibilities. The global χ 2 analysis using an extensive set of LHC data on charm meson R AA and v 2 has allowed us to constrain the preferred range of the two parameters: soft-hard intermediate cutoff |t * | = 1.5m 2 D and the scale of QCD coupling constant µ = πT . It is found that κ T /L have a mild sensitivity to t * , supporting the   Fig. 7 but for Pb-Pb collisions at √ sNN = 5.02 TeV. The data (solid [56], open [70]) are shown for comparison.
validity of the soft-hard approach when the coupling is not small. With this factorization model, we have calculated the transport coefficientq, drag coefficient η D , spatial diffusion coefficient 2πT D s , and then compared with other theoretical calculations and phenomenological extractions. Our analysis suggests that a small value 2πT D s 6 appears to be much preferred near T c . Finally we've demonstrated a simultaneous description of charm meson R AA and v 2 observables in the range p T ≤ 8 GeV. We've further made predictions for bottom meson observables in the same model.
We end with discussions on a few important caveats in the present study that call for future studies:  • In the framework of Langevin approach, the heavy flavor dynamics is encoded into three coefficient, κ T , κ L and η D , satisfying Einstein's relationship. It means that two of them are independent while the third one can be obtained accordingly. The final results therefore depend on the arbitrary choice of which of the two coefficients are calculated with the employed model [15,38]. For consistency, we cal-culate independently the momentum diffusion coefficients κ T /L with the factorization approach, and obtain the drag coefficient via η D = η D (κ T , κ L ) (see Eq. 4). A systematic study among the different options may help remedy this situation.
• According to the present modeling, the resulting spatial diffusion coefficient exhibits a relatively strong increase of temperature comparing with lattice QCD calculations, in particular at large T , as shown in Tab  • It is realized [72] that the heavy quark hadrochemistry, the abundance of various heavy flavor hadrons, provides special sensitivity to the heavylight coalescence mechanism and thus plays an important role to understand the observables like the baryon production and the baryon-to-meson ratio. A systematic comparison including the charmed baryons over a broad momentum region is therefore crucial for a better constraining of the model parameters, as well as a better extraction of the heavy quark transport coefficients in a model-todata approach.
• The elastic scattering (2 → 2) processes between heavy quark and QGP constituents are dominated for heavy quark with low to moderate transverse momentum [41]. Thus, here we consider only the elastic energy loss mechanisms to study the observables at p T ≤ 8 GeV. The missing radiative (2 ↔ 3) effects may help to reduce the discrepancy with R AA data in the vicinity of p T = 8 GeV, as mentioned above (see Fig. 8). With the soft-hard factorized approach, it would be interesting to include both elastic and radiative contributions in a simultaneous best fit to data in the whole p T region.
We also plan to explore this idea in the future.
In the QED case, µ + γ → µ + γ, one can calculate relevant interaction rate with small momentum transfer by using the imaginary part of the muon self-energy Σ(p 1 ) [73] Γ (A1) where, p 1 = (E 1 , p 1 ) and m 1 are the four-momentum and mass of the injected muon, respectively § . The trace term in Eq. A1 was calculated with a resummed photon propagator, which is very similar with the one in QCD ¶ . It is realized [74,75] that, in the QCD case, the contributions to Γ can be obtained from the corresponding QED calculations by simply substitution: e 2 → C F g 2 , where, e (g) is the QED (QCD) coupling constant and C F = 4/3 is the quark Casimir factor. It yields [74] Tr (p / 1 + m 1 )ImΣ(p 1 ) with where, v 1 = p 1 /E 1 denotes the HQ velocity; n B/F (E) = (e E/T ∓ 1) −1 indicates the thermal distributions for Bosons/Fermions andn B/F ≡ 1 ± n B/F accounts for the Bose-enhancement or Pauli-blocking effect. Here we have used the short notation q ≡ d 3 q (2π) 3 for phase space integrals. Taking | p 1 | and m 1/3 = m Q are both much greater than the underlying medium temperature, § The notations for the injected muon are same with the ones for the injected HQ in the elastic process. ¶ The structure of the QED HTL-propagator follows Eq. A10 but with opposite signs [48].
i.e. | p 1 |, m 1/3 T , thus, E 1/3 T and n F (E 3 ) is exponentially suppressed and can be dropped. Moreover, the first δ funciton in Eq. A4 can be simplified since Concerning the second δ funciton, it cannot contribute for ω less than or on the order of T [74], which will be deleted in this work. Finally, Eq. A3 and A4 can be reduced to By substituting Eq. A5a and A5b into Eq. A2, one gets can be rewritten as with the transverse and longitudinal spectral functions are given by the imaginary part of the retarded propagator We note that, in the weak coupling limit, a consistent method is to use the HTL resummed propagators, which is contributed by the quasiparticle poles and the Landau damping cuts [51]. In this analysis, we mainly focus on the low frequency excitation (|ω| < q), where the Landau damping is dominant and the quasiparticle excitation is irrelevant. The resulting spectral function is therefore denoted by ρ in Eq. A8.
The retarded propagator in Eq.A8 reads which is defined by setting q 0 = ω + iη (η → 0 + ), i.e. the real energy, for the dressed gluon propagator The medium effects are embedded in the HTL gluon selfenergy where, x = q 0 /q; Q(x) is the Legendre polynomial of second kind and m 2 D is the Debye screening mass squared for gluon The coupling constant, g, is quantified by the two-loop QCD beta-function [76] g −2 (µ) = 2β 0 ln( µ where, πT ≤ µ ≤ 3πT and Λ QCD = 261 MeV. N f is the number of active flavors in the QGP. Finally, for spacelike momentum, Eq. A8 can be expressed as with which Eq. A7 is computable. The momentum diffusion coefficients can be calculated by subsituting Eqs. A7, A15 and A16 back into Eq. 7 and 8, respectively, and then performing the angular integral, (A18) ¶ Note that the spectral functions (Eq. A15 and A16) are odd, and the resulting ρ T /L (−ω, q) = −ρ T /L (ω, q) are used to obtain Eq. A17 and A18.
with v 1 ≡ | v 1 |. The maximum momentum exchange is q max = √ 4E 1 T in the high-energy limit [77]. Next, we implement the further calculations by performing a simple change of variables, Using Eq. A19 and A20, one arrives at Eq. 9-12. Similar results can be found in Ref. [33,78].

Appendix B: Derivation of interaction rate and momentum diffusion coefficients in hard collisions
For two-boday scatterings, the transition rate is defined as the rate of collisions with medium parton i, which changes the momentum of the HQ (parton i) from p 1 ( p 2 ) to p 3 = p 1 − q ( p 4 = p 2 + q), ω Qi ( p 1 , q, T ) = p2 n(E 2 )n(E 3 )n(E 4 )v rel dσ Qi ( p 1 , p 2 → p 3 , p 4 ).
(B1) Typically, one can assumen(E 3 ) = 1 by neglecting the thermal effects on the HQ after scattering. The differential cross section summed over the spin/polarization and color of the final partons and averaged over those of incident partons, v rel dσ Qi ( p 1 , p 2 → p 3 , p 4 ) where, v rel = ( (p 1 · p 2 ) 2 − (m 1 m 2 ) 2 )/(E 1 E 2 ) is the relative velocity between the projectile HQ and the target parton. The interaction rate for a given elastic process reads Γ Qi (E 1 , T ) = d 3 q ω Qi ( p 1 , q, T ) (B1,B2) The momentum diffusion coefficients can be obtained by inserting Eq. B3 into Eq. 7 and 8, yielding (1) for hard collisions the momentum exchange is constrained by imposing θ(|t|−|t * |) in the first equality of Eq. B4 and B5; (2) the tree level matrix elements squared includes the contributions from the various channels, which are given in Ref. [79]: • for Q + q → Q + q (t-channel only) |M 2 | Qq (s, t) = 2N f · 2N c · g 4 4 9 (m 2 1 − u) 2 + (s − m 2 1 ) 2 + 2m 2 1 t t 2 (B6) • for Q + g → Q + g (t, s and u-channel combined) Concerning the degeneracy factors, in Eq. B6, 2N f reflects the identical contribution from all light quark and anti-quark flavors, and 2N c indicates the summing, rather than averaging, over the helicities and colors of the incident light quark, while in Eq. B7, the factor 2(N 2 c − 1) denotes the summing over the polarization and colors of the incident gluon. The running coupling constant takes g(µ) 2 = 4πα s (µ), which is given by Eq. A14a with the scale µ = √ −t.
(3) with the help of δ-function, we can reduce the integral in Eq. B4 and B5 from 9-dimension (9D) to 4D in the numerical calculations, by transforming the integration variables from ( p 2 , p 3 , p 4 ) to (| p 2 |, cosψ, t), where ψ is the polar angle of p 2 ; it yields the resutls as shown in the second equality of Eq. B4 and B5; the relevant limits of integration together with the additional notations are summarized below: (4) then, one can follow the procedure of Ref. [48] for the analytical evaluation of ω integral. The obtained results for κ T and κ L are shown in Eq. 13 and 14, respectively.