Parametric estimate of the relative photon yields from the Glasma and the Quark-Gluon Plasma in heavy-ion collisions

TL;DR

This paper analyzes photon production during different stages of heavy-ion collisions, highlighting the significant contribution of the Glasma phase relative to the thermalized Quark-Gluon Plasma, especially at high saturation scales.

Contribution

It provides a parametric study quantifying the photon yields from the Glasma and QGP stages, emphasizing the importance of non-equilibrium effects in photon production.

Findings

Glasma contributes significantly to photon yields at high saturation scales.

Non-equilibrium effects are crucial for accurate photon yield calculations.

Results align with experimental charged hadron multiplicities at RHIC and LHC.

Abstract

Recent classical-statistical numerical simulations have established the "bottom-up" thermalization scenario of Baier et al. as the correct weak coupling effective theory for thermalization in ultrarelativistic heavy-ion collisions. We perform a parametric study of photon production in the various stages of this bottom-up framework to ascertain the relative contribution of the off-equilibrium "Glasma" relative to that of a thermalized Quark-Gluon Plasma. Taking into account the constraints imposed by the measured charged hadron multiplicities at RHIC and the LHC, we find that Glasma contributions are important especially for large values of the saturation scale at both energies. These non-equilibrium effects should therefore be taken into account in studies where weak coupling methods are employed to compute photon yields.

Equations

$$E\frac{dN}{d^{4}X{d}^{3}p}=\frac{1}{2(2\pi {)}^3}{\int }_{p_1,p_2,p_3}|\mathcal{M}{|}^2(2\pi {)}^4{\delta }^4(P_1+P_2-P_3-P)f_1(p_1)f_2(p_2)\left[1\pm f_3(p_3)\right],$$

$${\int }_p=\int \frac{d^{3}p}{(2\pi {)}^{3}2E_p}$$

$$|{\mathcal{M}}_{\text{anni}}{|}^2=\frac{160}{9}16{\pi }^2\alpha {\alpha }_s\frac{u^2+t^2}{ut},$$

$$|{\mathcal{M}}_{\text{Comp}}{|}^2=\frac{320}{9}16{\pi }^2\alpha {\alpha }_s\frac{u^2+s^2}{-us},$$

$$f_1(E_1)f_2(E_2)\simeq e^{-(E_1+E_2)/T}.$$

$$E\frac{d{N}^{\text{th}}}{d^{4}x{d}^{3}p}=K\frac{5}{9}\frac{\alpha {\alpha }_s}{2{\pi }^2}{T}^2{e}^{-E/T}\log \left(\frac{2.912}{g^2}\frac{E}{T}\right).$$

$$d^{4}X=\tau d\tau d\eta d^2{x}_{\perp },$$

$$\frac{d^{3}p}{E}=d{y}_p{d}^2{p}_{\perp },$$

$$\frac{1}{S_{\perp }}\frac{d{N}^{\text{th}}}{d{y}_p{d}^2{p}_{\perp }}$$

$$\simeq \frac{5}{9}C\frac{\alpha {\alpha }_s}{2{\pi }^2}\int \tau d\tau T^2\int d\eta e^{-E/T},$$

$$C=K\log \left(1+\frac{2.912}{g^2}\right).$$

$$u_{\mu }=(\cosh \eta ,0,0,\sinh \eta ).$$

$${\int }_{-\infty }^{\infty }d\eta e^{-p^{\mu }{u}_{\mu }/T}$$

$$\frac{1}{S_{\perp }}\frac{d{N}^{\text{th}}}{d{y}_p{d}^2{p}_{\perp }}=\frac{10}{9}C\frac{\alpha {\alpha }_s}{2{\pi }^2}{\int }_{{\tau }_{\text{th}}}^{{\tau }_c}\tau d\tau T^2{K}_0\left(p_{\perp }/T\right).$$

$${\int }_0^{\infty }dxx{K}_0(x)=1,$$

$$\frac{1}{S_{\perp }}\frac{d{N}^{\text{th}}}{d{y}_p}=\frac{10}{9}C\frac{\alpha {\alpha }_s}{2{\pi }^2}{\int }_{{\tau }_{\text{th}}}^{{\tau }_c}\tau d\tau T^4.$$

$${\int }_1^{\infty }dxx{K}_0(x)=K_1(1)=\mathrm{0.601...}.$$

$${\partial }_{\tau }\mathcal{E}=-\frac{\mathcal{E}+P_L}{\tau }$$

$$P_L=\frac{1}{3}\mathcal{E},$$

$$\mathcal{E}(\tau )=\mathcal{E}({\tau }_{\text{th}}){\left(\frac{{\tau }_{\text{th}}}{\tau }\right)}^{4/3},$$

$$T(\tau )=T_{\text{th}}{\left(\frac{{\tau }_{\text{th}}}{\tau }\right)}^{1/3},$$

$$\frac{1}{S_{\perp }}\frac{d{N}^{\text{th}}}{d{y}_p}=\frac{5}{3}C\frac{\alpha {\alpha }_s}{2{\pi }^2}{T}_{\text{th}}^4{\tau }_{\text{th}}^2\left[{\left(\frac{{\tau }_c}{{\tau }_{\text{th}}}\right)}^{2/3}-1\right].$$

$${\tau }_{\text{th}}\simeq c_{\text{eq}}{\alpha }_s^{-13/5}{Q}_s^{-1}\text{and}{T}_{\text{th}}\simeq c_T{c}_{\text{eq}}{\alpha }_s^{2/5}{Q}_s.$$

$$\frac{1}{Q_s^2{S}_{\perp }}\frac{d{N}^{\text{th}}}{d{y}_p}\simeq \frac{5}{3}{c}_{\text{eq}}^6{c}_T^{4}C\frac{\alpha }{2{\pi }^2}{\alpha }_s^{-13/5}\left[{\left(\frac{{\tau }_c}{{\tau }_{\text{th}}}\right)}^{2/3}-1\right],$$

$$\frac{{\tau }_c}{{\tau }_{\text{th}}}={\left(\frac{T_{\text{th}}}{T_c}\right)}^3,$$

$$\frac{1}{Q_s^2{S}_{\perp }}\frac{d{N}^{\text{th}}}{d{y}_p}\simeq \frac{5}{3}{c}_{\text{eq}}^6{c}_T^{4}C\frac{\alpha }{2{\pi }^2}{\alpha }_s^{-13/5}\left[c_{\text{eq}}^2{c}_T^2{\alpha }_s^{4/5}{\left(\frac{Q_s}{T_c}\right)}^2-1\right].$$

$$E\frac{dN}{d^{4}X{d}^{3}p}=\frac{40}{9{\pi }^2}\alpha {\alpha }_s\mathcal{L}{f}_q(\mathbf{p})\int \frac{d^3{p}^{\prime }}{(2\pi {)}^3}\frac{1}{p^{\prime }}\left[f_g({\mathbf{p}}^{\prime })+f_q({\mathbf{p}}^{\prime })\right],$$

$$\mathcal{L}=\int \frac{dq}{q},,$$

$$\frac{1}{S_{\perp }}\frac{dN}{d{y}_p}=\frac{40}{9{\pi }^2}(2\pi {)}^3\alpha {\alpha }_s\mathcal{L}\int \tau d\tau I_q(\tau )\left[I_g(\tau )+I_q(\tau )\right],$$

$$I_{g/q}(\tau )=\int \frac{d^{3}p}{(2\pi {)}^3}\frac{1}{p}{f}_{g/q}(\mathbf{p},\tau ).$$

$$I_g(\tau )$$

Full text

I Introduction

Significant theoretical progress in understanding the space-time evolution of ultrarelativistic heavy-ion collisions can be achieved in the idealized high-energy limit where the QCD coupling $\alpha_{s}\ll 1$. The ab initio dynamics of such a system shortly after the collision corresponds to that of an over-occupied, strongly correlated non-Abelian plasma, exploding into the vacuum along the beam axis of the colliding nuclei Kovner:1995ja ; Krasnitz:1998ns ; Gelis:2007kn . The properties of this over-occupied plasma of strongly correlated quarks and gluons, often called a Glasma Lappi:2006fp , can be determined by employing classical-statistical methods.

An unfortunate complication is that the ab initio classical-statistical framework breaks down when the gluon occupancy becomes of order unity and the “quantum one-half” contributions become comparable to the leading classical contributions in real-time correlation functions Aarts:1997kp ; Berges:2013lsa . When this occurs, quantum kinetic descriptions are appropriate. However, because there is a significant overlap between the classical and quantum regimes in the evolution of the Glasma Mueller:2002gd ; Jeon:2004dh , classical-statistical simulations can help identify the right effective kinetic theory for the subsequent evolution of the Glasma to thermal equilibrium. The proper matching of the two frameworks is essential because the complexity of the dynamics of infrared modes in the system can lead to a number of weak coupling kinetic thermalization scenarios Baier:2000sb ; Blaizot:2011xf ; Kurkela:2011ti ; Bodeker:2005nv .

A recent breakthrough was achieved through large-scale numerical simulations of expanding non-Abelian plasmas in weak coupling where it was demonstrated that the Glasma flows to a non-thermal fixed point Berges:2013eia ; Berges:2013fga that is insensitive to details of the initial conditions111It is worth noting that Color Glass Condensate (CGC) Gelis:2010nm initial conditions for the Glasma Epelbaum:2013waa lead very rapidly to a gluon number over-occupancy Berges:2014yta that subsequently flows to this non-thermal fixed point.. Remarkably, the non-thermal fixed point identified by the classical-statistical simulations corresponds to the early stage of the “bottom-up” thermalization scenario of Baier et al. Baier:2000sb –henceforth referred to by the acronym BMSS. This result was unanticipated because it was believed previously Arnold:2003rq that plasma instabilities (not included in the BMSS framework) should in principle play a big role in kinetic realizations of the expanding Glasma Mrowczynski:1993qm ; Rebhan:2005re ; Romatschke:2005pm ; Mrowczynski:2016etf . The puzzling absence of late time plasma instabilities in the Glasma is strongly indicative of the large role of infrared modes as suggested by numerical results on the longitudinal to transverse pressure ratio Berges:2015ixa and in the striking universality of the non-Abelian non-thermal fixed point to that exhibited by expanding self-interacting scalar fields Berges:2014bba . Thus the BMSS kinetic theory–in the regime where occupancies are large–is best regarded as an effective description that captures the correct physics of the Glasma, in analogy to effective kinetic descriptions of weak wave turbulence Zakharov-book .

Given the significant developments we outlined in the context of real-time studies of early times, it is important to understand their phenomenological consequences for heavy-ion collisions. These can be ascertained by extrapolating the weak coupling results to realistic computations. For non-thermal fixed points in scalar theories, it has been shown recently that such extrapolations of the classical-statistical results are robust even for values of the scalar coupling constant $\lambda\sim 1$ Berges:2016nru . For gauge theories, the validity of such extrapolations is open to question; it is nevertheless useful to perform such an extrapolation and understand the phenomenological consequences thereof. Adopting this point of view, it has been shown recently Kurkela:2015qoa ; Keegan:2016cpi that a sophisticated implementation of the effective kinetic theory Arnold:2002zm can be smoothly matched to relativistic viscous hydrodynamics on the early time scales required by heavy-ion phenomenology Gale:2013da .

Electromagnetic signatures such as photon production are uniquely sensitive to the different stages of evolution in heavy-ion collisions. The relative rates of their production from the Quark-Gluon Plasma (QGP) and hadron gas stages of the evolution have been discussed for some time Kapusta:1991qp ; Shen:2013vja . The measurement of the elliptic flow coefficient for photons (which arise of course from the underlying anistropic flow of quarks and hadrons Chatterjee:2005de ) adds another handle to probe the space-time dynamics at the time of emission Adare:2015lcd ; Adam:2015lda . Comparisons of hydrodynamic and transport models, which implement contributions to the photon yields and elliptic flow from both thermal QGP emission and from the hadron gas, to the available data, indicate that reproducing both the photon yield and elliptic flow simultaneously is challenging in these models Chatterjee:2013naa ; Shen:2014lpa ; vanHees:2014ida ; Bratkovskaya:2014mva ; Paquet:2015lta . This discrepancy has been dubbed the “direct photon puzzle”–for a recent review, see Shen:2016odt .

The above estimates do not include the contribution from pre-equilibrium photon production222For a discussion in strong coupling frameworks, see Iatrakis:2016ugz . For a recent discussion in the context of partonic transport simulations, see Greif:2016jeb .. These may be especially important in semi-peripheral heavy-ion collisions and in proton-nucleus collisions, where their relative contribution to the thermal yield Shen:2015qba , as a function of centrality, may help constrain the onset of thermalization in QCD matter. First phenomenological studies of photon production in the Glasma McLerran:2014hza ; McLerran:2014oea suggest that this contribution is significant. In particular, it is argued McLerran:2015mda that photon yields and elliptic flow can be explained if the Glasma has a hard component which thermalizes at late times relative to thermalization times typically assumed in hydrodynamic simulations.

In this work, we will estimate the rate of photon production within the BMSS framework. In particular, we will extract the parametric dependence of the rate, from different stages of the evolution of the Glasma, on the QCD coupling constant. These rates will be compared to that of photon production from the equilibrated QGP. We find that our conclusions are very sensitive to the saturation scale $Q_{s}$ in the Glasma and coefficients that relate this scale to the initial temperature and thermalization time. These coefficients can be fixed to fair accuracy from heavy-ion data on charged particle multiplicites; it is nevertheless subject to significant systematic uncertainties that we will elaborate on.

As noted, the results quoted in this paper are parametric estimates and valid only in the kinetic regime. The results are however interesting enough to suggest that a more detailed computation of pre-equilibrium photon production rates is desirable. When gluon occupancies are large, this can be achieved by using classical-statistical simulations including dynamical quarks Gelis:2005pb ; Gelis:2015eua ; Gelfand:2016prm 333Similar computations have been performed in proton-nucleus collisions Gelis:2002ki ; Benic:2016yqt ; Benic:2016uku . In this case, at least for minimum-bias collisions, Glasma evolution is not significant.. First results for the photon production from color background fields in a fixed box geometry are now available Tanji:2015ata . When the occupancy of gluons becomes of order unity and smaller, one needs to solve the coupled set of Boltzmann equations for both quark and gluon evolution, with the initial conditions given by the results of the classical-statistical simulations. Again, for a fixed box geometry, this computation has been done Blaizot:2014jna . In an accompanying paper, we will discuss the extension of these results to that of the expanding Glasma Tanji:2017suk .

The paper is organized as follows. In the next section, we will outline the well known result for photon production from a thermal QGP. Unlike previous computations though, we will study the parametric dependence of the integrated rate on the QCD coupling taking into account the fact that the thermalization time too depends on the coupling. The corresponding estimates for the various stages of Glasma evolution are given in section III. In section IV, we will first discuss how one can constrain the numerical coefficients in the initial temperature and thermalization time from data on charged particle multiplicities. Given these constrained values, we will then compute the relative yields for photon production from the Glasma and the QGP for varying system size and center-of-mass energy. We will end with a summary and an outlook on future work. A formula for the photon production rate in the small-angle approximation is derived in Appendix A.

II The thermal photon production rate

We will begin this section with a brief recapitulation of well known computations in the literature on the thermal photon production rate. We will subsequently embed this rate in an expanding geometry and discuss the consequences of the parametric dependence of the integrated rate on the BMSS thermalization time relative to the hadronization time.

II.1 Recap of known thermal production estimate

We will rely on the kinetic expression for the production of on-shell photons with momentum ${\boldsymbol{p}}=(p_{x},p_{y},p_{z})$ at the space-time point $X=(t,x,y,z)$ from two-to-two scattering Kapusta:1991qp ; Baier:1991em ,

[TABLE]

where

[TABLE]

and $P=(E_{p},{\boldsymbol{p}})$. The squared amplitude $|\mathcal{M}|^{2}$ should be understood as summed over spins, colors and flavors of all incoming and outgoing particles. For massless up and down quarks, the amplitude is given in terms of the Mandelstam variables $s=(p_{1}+p_{2})^{2}$, $t=(p_{1}-p)^{2}$ and $u=(p_{3}-p_{1})^{2}$ as

[TABLE]

for the annihilation process, and

[TABLE]

for Compton scattering with the electromagnetic coupling constant $\alpha$. Because photons are never in equilibrium in the QGP, the formula (1) is valid for both equilibrium and non-equilibrium rates as long as the kinetic description is applicable.

For photon production from a thermal medium, we can simplify subsequent computations considerably by adopting the approximation

[TABLE]

This is a good approximation for high-energy photons with $E_{1}+E_{2}>E\gg T$. For $f_{3}(E_{3})$, one has to keep the Bose-Einstein or the Fermi-Dirac distribution form since $E_{3}$ is not necessarily large. Using this simplification, one can arrive at the formula for the thermal photon production rate from two-to-two scattering Kapusta:1991qp ,

[TABLE]

In the derivation of this expression, the infrared divergence is regulated by the Hard Thermal Loop resummed quark propagator. We have introduced a constant $K$, which we will explain in the next paragraph.

This rate is leading order (LO) in $\alpha$ and $\alpha_{s}$. However the two-to-three bremsstrahlung and pair annihilation processes, which are naively higher order, contribute at parametrically the same order as the two-to-two processes Aurenche:1998nw ; Aurenche:1999tq ; Aurenche:2000gf . The complete LO results444We note that the next-to-leading $\mathcal{O}(g)$ contribution to the thermal photon rate has been computed; this gives a surprisingly small correction to the LO rate even for $\alpha_{s}=0.3$ Ghiglieri:2013gia . which include the Landau-Pomeranchuk-Migdal effect, as well as these collinear enhanced processes, have been derived by Arnold, Moore and Yaffe Arnold:2001ms . For $\alpha_{s}=0.2$, the naive LO rate differs from the complete LO rate by a factor of two in the photon momentum range $2.5\leq k/T\leq 10$. We introduced the $K$ factor in the above equation to approximately take these effects into account. For the purposes of this study, we will assume $K\simeq 2$ and shall employ the LO formula Eq. (6) henceforth to estimate the thermal photon yield.

II.2 BMSS estimate of the thermal production rate

With Eq. (6) in hand, we shall now make a parametric estimate of the photon yield in the thermal QGP stage based on a simple model for the space-time evolution and the BMSS results on the thermalization time and the initial temperature of the system Baier:2000sb .

The longitudinally expanding system is conveniently described by the comoving coordinates $\tau=\sqrt{t^{2}-z^{2}}$, $\eta=\text{arctanh}(z/t)$, ${\boldsymbol{x}}_{\!\perp}=(x,y)$ and momentum variables $y_{p}=\text{arctanh}(p_{z}/E)$, ${\boldsymbol{p}}_{\!\perp}=(p_{x},p_{y})$. We will assume for simplicity that the expansion is boost-invariant. The temperature is then a function of proper time, $T=T(\tau)$. By using the relations

[TABLE]

and

[TABLE]

we can rewrite (6) as

[TABLE]

with $S_{\perp}$ denoting the transverse area. In the second line, we have assumed that the $\tau$ and $\eta$ dependence of the logarithmic factor is slower than other factors, and replaced the logarithmic function555We have added 1 to the argument of the logarithm–as noted in Kapusta:1991qp , adding this constant gives very good agreement between the numerical computation of the rate and the analytic approximation from the $E\gg T$ region where it is justified, down to $E\simeq T$. $\log\left(1+\frac{2.912}{g^{2}}\frac{E}{T}\right)$ by a numerical factor $\log\left(1+\frac{2.912}{g^{2}}\right)$, which is denoted by $C$ after being combined with the $K$ factor:

[TABLE]

In a Lorentz covariant description, the photon energy $E$ should be replaced by $p^{\mu}u_{\mu}$, where $p^{\mu}$ is the energy-momentum 4-vector and $u_{\mu}$ is the co-moving 4-velocity,

[TABLE]

The $\eta$ integration can be done analytically,

[TABLE]

where $K_{n}(z)$ is the modified Bessel function of the second kind. To estimate the QGP photon yield, we need to define $\tau_{c}$, the time until which any weak coupling estimate of photon production in the QGP may at all be applicable. The most optimistic estimate of $\tau_{c}$ is the “hadronization time”, the time at which the temperature of the system reaches the crossover temperature, namely $T=T_{c}$. We shall adopt this definition henceforth. The photon yield in the QGP can then be expressed as

[TABLE]

By using the formula

[TABLE]

we can further integrate the yield over the transverse momentum,

[TABLE]

If we integrate $p_{\!\perp}$ from $T$ to $+\infty$ since the formula (6) is valid for p_{\!\perp}\raisebox{-2.15277pt}{>\stackrel{{\scriptstyle\textstyle>}}{{\sim}}>}T, we obtain a numerical factor of

[TABLE]

This uncertainty can be absorbed by the factor $C$, which anyway possesses order-one uncertainty.

In order to estimate the thermal photon yield (15), we need a profile for the time dependence of the temperature. For simplicity, we assume a 1+1 dimensional hydrodynamic expansion,

[TABLE]

with the equation of state

[TABLE]

where $\mathcal{E}$ is the energy density and $P_{L}$ is the longitudinal pressure. These equations can be easily solved to give

[TABLE]

and from the relation $\mathcal{E}\propto T^{4}$, the $\tau$ dependence of the temperature is found to be

[TABLE]

with $T_{\text{th}}=T(\tau_{\text{th}})$ being the temperature at the thermalization time. Plugging these expressions into Eq. (15) leads to the result

[TABLE]

We now observe that in BMSS Baier:2000sb ; Baier:2002bt , the time it takes the system to thermalize and the initial temperature are derived respectively to be

[TABLE]

Here $c_{\text{eq}}$ is a constant denoting the uncertainty in the BMSS estimate of the time it takes for the Glasma to thermalize – it can in principle be determined self-consistently in the BMSS framework. As we shall discuss, we will constrain it with data from RHIC and the LHC on hadron multiplicities. There is an additional constant $c_{T}$ which is needed to determine the initial temperature of the QGP. BMSS determine this number to be, to logarithmic accuracy, $0.16c$, where $c$ is the gluon liberation constant first discussed in Mueller:1999fp . The coefficient $c$ is a measure of how efficiently gluons in the wavefunction are released in the collision and can be estimated using boost-invariant classical Yang-Mills simulations of the Glasma Krasnitz:2000gz ; Krasnitz:2003jw ; Lappi:2003bi ; the most sophisticated estimate Lappi:2007ku gives this gluon liberation coefficient to be $c=1.1$. Thus $c_{T}\simeq 0.18$; however since there are additional logarithmic uncertainties, and $c$ is not known for the full 3+1-D Yang-Mills simulations, we will treat $c_{T}$ as a constant to be varied within a factor of two of the BMSS value.

Substituting these expressions for the thermalization time and the initial temperature in our expression for the rate, we obtain

[TABLE]

The ratio of the two time scales in the above expression can be obtained from the temporal profile of the temperature (20) as

[TABLE]

where $T_{c}=154\pm 9$ MeV is the crossover temperature in QCD Borsanyi:2010cj ; Bazavov:2011nk ; for our study, we will simply take $T_{c}=154$ MeV. Substituting this back in Eq. (23), we obtain our final expression for the thermal yield to be

[TABLE]

We now summarize the several undetermined constants in our weak coupling expression for the thermal photon yield. The constant $C$, given by Eq. (10), is simply the $K$-factor from uncertainities in the thermal rate modulo the logarithmic contribution. As discussed, we will henceforth take $K=2$. The constants $c_{\text{eq}}$ and $c_{T}$ cannot be completely determined from theory at present. As we will discuss in Sec. IV.1, these constants can be constrained using data on charged hadron multiplicities measured at RHIC and the LHC. Finally, for $\alpha_{s}$, we will take the one-loop value assuming it runs with the scale $Q_{s}$. (This choice also has systematic uncertainties which should be taken into account.) Modulo the stated uncertainties, our result for the thermal photon yield is a function of $Q_{s}$ alone, which varies both with system size and center-of-mass energy.

III Estimate of pre-equilibrium photon production

In the bottom-up Glasma thermalization scenario of BMSS Baier:2000sb , the pre-equilibrium Glasma evolution is divided into three temporal stages:

• (i)

$Q_{s}^{-1}\ll\tau\ll Q_{s}^{-1}\alpha_{s}^{-3/2}$

• (ii)

$Q_{s}^{-1}\alpha_{s}^{-3/2}\ll\tau\ll Q_{s}^{-1}\alpha_{s}^{-5/2}$

• (iii)

$Q_{s}^{-1}\alpha_{s}^{-5/2}\ll\tau\ll Q_{s}^{-1}\alpha_{s}^{-13/5}$ .

In stage (i), the gluons are highly occupied, ranging from an occupancy of $f\sim 1/\alpha_{s}$ at $\tau\sim 1/Q_{s}$ to unity at $\tau\sim Q_{s}^{-1}\alpha_{s}^{-3/2}$. The occupancy decreases in time as $(Q_{s}\tau)^{-2/3}$ as a consequence of the broadening of the longitudinal momentum distribution by elastic scatterings amongst hard gluons. These modify the typical longitudinal momentum from $p_{z}\sim 1/\tau$ to $p_{z}\sim\tau^{-1/3}$. This BMSS prediction is confirmed by the classical-statistical lattice simulations Berges:2013lsa ; Berges:2014yta which are a good approximation to the real-time dynamics of the theory for $f\gg 1$.

The quantum kinetic dynamics of the BMSS framework Arnold:2002zm underlies the dynamics of the stages (ii) and (iii), where the occupancy of hard gluons is less than unity. In stage (ii), the soft gluons that are being produced as a result of inelastic scattering dominate screening by providing a larger contribution to the Debye mass relative to that of hard gluons. In stage (ii), the typical longitudinal momentum of hard gluons is $p_{z}^{2}\sim\alpha_{s}\,Q_{s}^{2}$, which does not depend on time anymore. The anisotropy thus saturates at a value of the ratio of longitudinal to transverse pressure $P_{L}/P_{T}\sim\alpha_{s}$. The multiplicity of soft gluons is however still significantly smaller than those of hard gluons.

This is no longer the case by the start of stage (iii), with soft gluons dominating the multiplicity for $\tau>Q_{s}^{-1}\alpha_{s}^{-5/2}$. Further, $\tau>\tau_{\rm rel}$, the relaxation time of soft gluons; this indicates that the soft gluons have thermalized by then. The hard gluons however are not thermal, and they thermalize by losing energy to the heat bath of soft gluons through a process which corresponds to the description of jet quenching Baier:2000mf . The infusion of energy into the heat bath raises its temperature (even though the system is expanding) temporarily to saturate finally at the previously quoted temperature of $T_{\text{th}}=c_{T}\,c_{\text{eq}}\,\alpha_{s}^{2/5}Q_{s}$ and at the thermalization time of $\tau_{\text{th}}=c_{\text{eq}}\,\alpha_{s}^{-13/5}Q_{s}^{-1}$. Subsequently, the system undergoes hydrodynamical expansion, with the temperature of the system cooling as $T\sim\tau^{-1/3}$.

We will now estimate photon production from the three stages of the pre-equilibrium evolution of the Glasma.

III.1 Glasma stage (i)

To compute the scattering rate in Eq. (1), we will explicitly employ the small-angle approximation landau1936kinetic ; lifshitz1981physical , which dominates the $2\leftrightarrow 2$ scattering of energetic partons. By this approximation, the photon production rate is simplified to

[TABLE]

where $f_{g}$ and $f_{q}$ denote the momentum distribution functions of gluons and quarks, respectively. The details of this derivation are given in Appendix A. The symbol $\mathcal{L}$ denotes the so-called Coulomb logarithm,

[TABLE]

which should be regularized by infrared and ultraviolet cutoffs. In thermal equilibrium, these are respectively the Debye mass $m_{D}$ and the temperature, giving ${\cal L}\sim\log(1/g)$, which is the origin of the logarithm in Eq. (6).

By integrating (26) over ${\boldsymbol{p}}$, one obtains the photon yield per unit rapidity

[TABLE]

where we have introduced integrals

[TABLE]

Now the problem is reduced to the evaluation of the integrals $I_{g}$ and $I_{q}$, which are much simpler than the original multi-dimensional integral in Eq. (1). In the first stage of the Glasma evolution, $I_{q}$ in the square bracket is negligible compared to $I_{g}$.666This indicates that the pair annihilation process is negligible compared with Compton scattering in the gluon dominated medium. These integrals are related to the Debye screening mass $m_{D}$ by the expression $m_{D}^{2}=4\,g^{2}(N_{c}I_{g}+N_{f}I_{q})$.

In the first stage of bottom-up thermalization, the Debye mass is dominated by hard gluons whose typical transverse momentum is $p_{\!\perp}\simeq Q_{s}$. The integral $I_{g}$ therefore can be approximately related to the number density of hard gluons $n_{\text{hard}}$ as

[TABLE]

where $2\,(N_{c}^{2}-1)$ is the degeneracy factor for gluons. Since the total number of hard gluons is approximately conserved at this stage, the number density $n_{\text{hard}}$ decreases in time as $1/\tau$. Hence the integral $I_{g}$ can be represented as

[TABLE]

where $\kappa_{g}$ is a dimensionless constant, which we will fix later.

This functional form of $I_{g}(\tau)$ is consistent with the following scaling behavior for the gluon distribution, confirmed by classical-statistical simulations Berges:2013lsa ; Berges:2014yta :

[TABLE]

Plugging this scaling form into Eq. (29) leads to the expression,

[TABLE]

If we identify the product of the two dimensionless integrals in the last expression with $\kappa_{g}$, this expression is equivalent to (31).

Since quarks undergo the same scattering processes as gluons, namely small-angle elastic collisions, it is natural to assume that the quark occupation number for hard modes with $p_{\!\perp}\sim Q_{s}$ decreases with the same power law as that of gluons,

[TABLE]

This functional form of the quark distribution at the scale $p_{\!\perp}\sim Q_{s}$ has been confirmed by us in a kinetic treatment a la BMSS Tanji:2017suk . We note that while quark pairs are produced copiously at early times, the gluon fusion process will continue to produce pairs as the system evolves, albeit at a diminishing rate, with this contribution amenable to a perturbative treatment in a kinetic approach Blaizot:2014jna ; Ruggieri:2015tsa . In weak coupling, this effect is negligible and the evolution of the quark spectrum is well described by the scaling form (34) Tanji:2017suk . Performing the same computation as in the gluon case, we can express the integral for quarks $I_{q}$ as

[TABLE]

where another constant $\kappa_{q}$ has been introduced. We note that the factor $1/\alpha_{s}$ in $I_{g}$ is absent here.

Substituting the expressions (31) and (35) into Eq. (28), we can express the photon yield in stage (i) as

[TABLE]

We have used here $\tau_{0}\sim Q_{s}^{-1}$ and $\tau_{1}\sim Q_{s}^{-1}\alpha_{s}^{-3/2}$. This estimate has a systematic uncertainty since we do not fix the numerical coefficients of $\tau_{0}\sim Q_{s}^{-1}$ and $\tau_{1}\sim Q_{s}^{-1}\alpha_{s}^{-3/2}$. However, the uncertainty is small because $\tau_{0}$ and $\tau_{1}$ both only appear inside the logarithmic factor.

We therefore have a good estimate of the non-equilibrium photon yield in the first stage of the Glasma evolution. In weak coupling, this yield is parametrically much smaller than the thermal yield in Eq. (25). However as we noted previously, the thermal rate is sensitive to high powers of the constants $c_{T}$ and $c_{\text{eq}}$. We will demonstrate later that, for realistic values of the coupling, the ratio of the two yields will depend strongly on these coefficients.

The values of the constants $\kappa_{g}$ and $\kappa_{q}$, which depend on the normalization of the scaling function $f_{S}$, can be determined by the following considerations. The number density of hard gluons produced immediately after the heavy-ion collision can be expressed as Mueller:1999fp

[TABLE]

where $c$ is the gluon liberation coefficient we discussed previously, with $c=1.1$ from solutions of the boost-invariant classical Yang-Mills equations Lappi:2007ku . Combining this with Eq. (30) immediately leads to

[TABLE]

To compute $\kappa_{q}$, we assume that the quark number density at $\tau\sim Q_{s}^{-1}$ is smaller than that of gluons by the factor $\alpha_{s}$. Multiplying $\alpha_{s}$ as well as a factor to convert the degeneracy factor to the gluon number density in Eq. (37), we obtain

[TABLE]

from which we find

[TABLE]

III.2 Glasma stage (ii)

In this stage, the typical occupancy of hard gluons is less than unity and decreases as $f\sim\alpha_{s}^{-3/2}/(Q_{s}\tau)$. While the total particle number is still dominated by hard gluons, the Debye mass is dominated by soft gluons, and it behaves as $m_{D}\sim\alpha_{s}^{3/8}Q_{s}(Q_{s}\tau)^{-1/4}$ Baier:2000sb . Therefore, the time dependence of the integral $I_{g}$ is different from its dependence in the Glasma first stage, and it can be evaluated as

[TABLE]

By assuming this expression agrees with that in the first stage777Since we have considered only the hard contribution to $I_{g}$ in stage (i) and the soft one in stage (ii), the two expressions for $I_{g}$ do not necessarily agree at $\tau_{1}$. However this assumption may be appropriate because, as shown in Baier:2000sb , the Debye mass receives equal contributions from hard and soft gluons at $\tau_{1}$., given by Eq. (31), at $\tau_{1}=c_{1}\,Q_{s}^{-1}\alpha_{s}^{-3/2}$, the overall normalization of $I_{g}$ is identified to be

[TABLE]

Here we have introduced the numerical coefficient $c_{1}$ to fix the time scale $\tau\sim Q_{s}^{-1}\alpha_{s}^{3/2}$.

For quarks, soft modes never dominate the integral $I_{q}$. Therefore we assume that the time dependence of $I_{q}$ is the same as that in the first stage888 This assumption may not be adequate. For chemical equilibration between quarks and gluons, the quark production process is essential. If the total quark number is increasing, the time dependence of $I_{q}$ should be slower than $\tau^{-1}$ and the photon yield may become larger than the present estimate.:

[TABLE]

Plugging these expressions into (28), we obtain

[TABLE]

To obtain the last expression, we have substituted $\tau_{1}=c_{1}\,Q_{s}^{-1}\alpha_{s}^{-3/2}$ and $\tau_{2}=c_{2}\,Q_{s}^{-1}\alpha_{s}^{-5/2}$, where $c_{2}$ is another numerical coefficient. Although the numerical coefficients $c_{1}$ and $c_{2}$ are unknown, it is reasonable to assume that they are both $\mathcal{O}(1)$ and their ratio is close to one. Because of the exponent 1/2, the photon yield is less sensitive to the uncertainty of these coefficients compared with the case of the thermal photon yield (25), in which the numerical coefficients $c_{\text{eq}}$ and $c_{T}$ appear in high powers. Here we simply assume that $c_{1}=c_{2}$. These coefficients can be determined more accurately by real-time lattice simulations and more detailed kinetic theory computations, which are beyond the scope of the present work. Finally, we obtain

[TABLE]

In the weak coupling limit, this yield is much larger than the photon yield from the first stage given in Eq. (36).

III.3 Glasma stage (iii)

In this final stage of Glasma evolution, the total particle number is dominated by soft gluons, while most of the energy is carried by a few hard particles. In the bottom-up thermalization scenario, the soft gluons reach thermal equilibrium first and the hard gluons subsequently lose their energy to the heat bath of soft gluons by a turbulent bremsstrahlung process Blaizot:2012fh . The temperature of the heat bath increases in time as

[TABLE]

Photons produced from this thermal bath can be estimated from the thermal rate (15) by replacing $\tau_{\text{th}}$ by $\tau_{2}=c_{2}\,Q_{s}^{-1}\alpha_{s}^{-5/2}$ and $\tau_{c}$ by $\tau_{\text{th}}=c_{\text{eq}}\,Q_{s}^{-1}\alpha_{s}^{-13/5}$:

[TABLE]

The second term inside the brackets has a systematic uncertainty since we do not fix the numerical coefficient $c_{2}$. However, since the temperature increases in this stage of Glasma evolution, the photon yield is dominated by $\tau\sim\tau_{\text{th}}$. Further, in the limit of weak coupling, the second term is negligible compared to the first term. In keeping with our assumption about coefficients multiplying time scales being of the same order, we shall simply henceforth replace the ratio $c_{2}/c_{\text{eq}}$ by unity.

It is interesting to compare the above expression for the Glasma yield in stage (iii) directly with the thermal yield in Eq. (25). We obtain

[TABLE]

In the limit of very weak coupling, which also corresponds to very large values of $Q_{s}$ (taking into account the running of the coupling with $Q_{s}$), the ratio of the two yields goes to zero. We will discuss this ratio for realistic values of $Q_{s}$, $c_{T}$ and $c_{\text{eq}}$ in Sec. IV.3.

One can also have photon emission during the process of quenching the hard quarks to the heat bath. This has been computed previously for a static medium Zakharov:2004bi for the case where only photons are radiated; since the contribution of soft gluon radiation along with that of photons may be significant, this analysis is incomplete. We leave further discussion of this contribution to future work.

IV Phenomenological estimates: Thermal versus Glasma photon yields

We will now obtain estimates for the thermal photon yield and the corresponding Glasma photon yield, based on the expressions in the previous sections. The biggest uncertainties are the parameters $c_{\text{eq}}$ and $c_{T}$ which appear, with high powers, in the thermal photon yield. As we shall now discuss, they are constrained by data from RHIC and the LHC on charged particle multiplicity.

IV.1 Estimates of $c_{\rm eq}$ and $c_{T}$

In Sec. II, we observed that the thermal photon yield is sensitive to the coefficients $c_{\rm eq}$ and $c_{T}$ that appear in estimates of the thermalization time and the initial temperature at that time. Assuming that the system satisfies nearly ideal hydrodynamic flow conserving entropy subsequently999It has been shown in Ref. Kurkela:2015qoa that weak-coupling dynamics matches smoothly with hydrodynamic simulations for the couplings we employ in our work. The viscous hydrodynamic simulations employing a temperature-dependent $\eta/s$ produce entropy, however, it is an effect of about 15% Schenke:2013dpa ; McDonald:2016vlt , which is part of the systematic uncertainties of our parametric estimates., one can use thermodynamic relations and the data on charged particle multiplicities to constrain them Baier:2002bt ; Baier:2011ap .

The entropy of hadrons per unit rapidity can be related to the measured multiplicity of charged hadrons as

[TABLE]

The proportionality constant $k_{S/N}$ can be phenomenologically estimated in several ways as summarized in Gubser:2008pc . We shall adopt the value $k_{S/N}=7.2$ that has been extracted from experimental data for particle yields, spectra and source sizes estimated by two-particle interferometry Pal:2003rz . The entropy of the QGP per unit rapidity at time $\tau$ and temperature $T$ is

[TABLE]

where $\nu_{\text{QGP}}=2\,(N_{c}^{2}-1)+\frac{7}{2}\,N_{f}\,N_{c}$ denotes the internal degrees of freedom for the QGP phase; $\nu_{\rm QGP}=37$ for $N_{f}=2$.

Equating the entropy in the two phases gives

[TABLE]

Since entropy is conserved in ideal hydrodynamic evolution, we can estimate the left hand side at any time during this stage. By substituting $\tau=\tau_{\text{th}}$ and $T=T_{\text{th}}$ given by Eq. (22), we can relate the unknown constants $c_{\text{eq}}$ and $c_{T}$ to the charged hadron multiplicity as

[TABLE]

We have introduced the number of participants $N_{\text{part}}$ because experimental data for $\frac{2}{N_{\text{part}}}\frac{dN_{\text{ch}}}{d\eta}$ are available as a function of $N_{\text{part}}$ at RHIC and the LHC Adler:2004zn ; Abelev:2008ab ; Aamodt:2010cz ; Alver:2010ck . In this section, we will express all quantities as a function of $N_{\text{part}}$. For that, we need a model which can relate $Q_{s}$ and $S_{\perp}$ to $N_{\text{part}}$.

In the IP-Glasma model Schenke:2012wb ; Schenke:2012fw ; Schenke:2013dpa , the IP-Sat dipole framework of gluon saturation Kowalski:2003hm ; Rezaeian:2012ji and the geometrical cross-sections of the Glauber model Miller:2007ri are combined, and $Q_{s}^{2}\,S_{\perp}\equiv\int dx_{T}^{2}\,Q_{s}^{2}(x_{T},\sqrt{s})$ can be determined as a function of $N_{\text{part}}$ and the center-of-mass energy $\sqrt{s}$.101010We thank Prithwish Tribedy for providing us with the values of $Q_{s}^{2}S_{\perp}$ and $N_{\text{part}}$ in the IP-Glasma model. We further use the Glauber model to compute $S_{\perp}$ as a function of $N_{\text{part}}$, and can compute $Q_{s}^{2}$ by combining the two results. However in the IP-Glasma model, there are two spatially varying saturation scales, one from the projectile and other from the target. At any given spatial position, it is the lower of the two saturation scales that governs the multiplicity of gluons produced locally. Thus the value of $Q_{s}^{2}$ extracted in the aforementioned manner is lower than the value that governs the typical momentum of the produced gluons. Thus while we will adopt the $N_{\text{part}}$ dependence of $Q_{s}$ provided by the IP-Glasma model, we will account for the harder momentum distribution of produced gluons by multiplying the IP-Glasma $Q_{s}$ by a numerical factor.

In Fig. 1, we plot the values of $Q_{s}^{2}$ for $\sqrt{s_{\text{NN}}}=200$ GeV at RHIC and $\sqrt{s_{\text{NN}}}=2.76$ TeV at the LHC. We have adjusted the overall normalization to obtain a “reference” $Q_{s}^{2}$ at the RHIC most central collision ($N_{\text{part}}=353$) of 2 GeV2. Absent a first principles determination of the hard scale in the Glasma, this choice of this reference value is somewhat arbitrary with the only consideration being that it is a semi-hard scale. The LHC values are obtained by multiplying the factor $(2.76/0.2)^{0.3}$ to the RHIC values at the same $N_{\text{part}}$. This energy dependence is consistent with the empirically observed value for particle multiplicites from RHIC to LHC energies. Because our simple estimate of photon production is less reliable for peripheral collisions, we do not plot the region $N_{\text{part}}<150$. In Sec. IV.3, we shall discuss the dependence of the photon yields for the stated $Q_{s}^{2}$ reference scale.

For the running coupling, we employ the one-loop expression,

[TABLE]

where $\Lambda_{\text{QCD}}=0.2$ GeV and we have assumed $N_{f}=3$. For the values of $Q_{s}^{2}$ plotted in Fig. 1, the variation of the coupling is rather slow; it varies between $\alpha_{s}=0.35$–0.42 for the $N_{\text{part}}$ range for $\sqrt{s_{\text{NN}}}=200$ GeV at RHIC and $\alpha_{s}=0.29$–0.33 for the comparable range for $\sqrt{s_{\text{NN}}}=2.76$ TeV at the LHC.

For given $Q_{s}$ and $S_{\perp}$, we can compute $c_{\text{eq}}\,c_{T}^{3/4}$ from the constraint relation in Eq. (52). For the values of charged hadron multiplicity, we have used the PHENIX data Adler:2004zn and the ALICE data Aamodt:2010cz . The results are shown in Fig. 2. The values of $c_{\text{eq}}\,c_{T}^{3/4}$ are almost independent of $N_{\text{part}}$ as they should be for our analysis to be self consistent. For $c_{T}=0.18$, the value of $c_{\text{eq}}$ ranges from 1.0 to 1.3.

Since we have only one constraint equation given by Eq. (52), the two coefficients $c_{\text{eq}}$ and $c_{T}$ cannot be fixed independently. Only the combination $c_{\text{eq}}\,c_{T}^{3/4}$ is fixed111111One may hope to fix the two coefficients independently by further using the measured transverse energy. However, the constraint from the transverse energy density is not independent of the hadron multiplicity constraint when the speed of sound is $c_{s}=1/\sqrt{3}$. If $c_{s}\neq 1/\sqrt{3}$, the two constraints are independent and the two coefficients can be fixed individually. However, the result is very sensitive to the value of $c_{s}$ and therefore involves a large uncertainty.. In the BMSS papers Baier:2002bt ; Baier:2011ap , $c_{T}\simeq 0.18$ is estimated121212Note that we have incorporated the gluon liberation constant $c\simeq 1.1$ in our estimate. to logarithmic accuracy. To indicate the impact of the uncertainty in this quantity, we will vary it by a factor of two in the range $c_{T}=0.1$–0.4.

IV.2 Estimates of $\tau_{\rm th}$, $T_{\rm th}$, and $\tau_{c}$ in the bottom-up thermalization scenario

Before the discussion on the photon yields, it is instructive to show the estimation of the thermalization time $\tau_{\text{th}}$, the initial temperature $T_{\text{th}}$, and the hadronization time $\tau_{c}$ as a function of $N_{\text{part}}$. Since the coefficients $c_{\text{eq}}$ and $c_{T}$ are constrained by the observed hadron multiplicities, we can numerically evaluate these quantities and compare those for RHIC and LHC energies. All the results shown in this subsection assume the $Q_{s}^{2}$ profile plotted in Fig. 1.

In Fig. 3, the thermalization time $\tau_{\text{th}}=c_{\text{eq}}\,\alpha_{s}^{-13/5}Q_{s}^{-1}$ is plotted as a function of $N_{\text{part}}$. One may expect that $\tau_{\text{th}}$ is a decreasing function of $N_{\text{part}}$ because $Q_{s}$ increases for increasing $N_{\text{part}}$. However, it is not always true because there is a competition between the factors $\alpha_{s}^{-13/5}$ and $Q_{s}^{-1}$. The running coupling $\alpha_{s}$ is a decreasing function of $N_{\text{part}}$. Although the variation of the coupling is slow, the factor $\alpha_{s}^{-13/5}$ varies relatively strongly and it tends to cancel the variation of $Q_{s}^{-1}$. This is the reason why the plots of $\tau_{\text{th}}$ are rather flat as a function of $N_{\text{part}}$ and the values of $\tau_{\text{th}}$ are similar for RHIC and LHC energies. By the same reason, $\tau_{\text{th}}$ is rather insensitive to the choice of the normalization for the $Q_{s}^{2}$ profile in a realistic parameter range. In the figure, the uncertainty of $\tau_{\text{th}}$ corresponding to the variation of $c_{T}=0.1$–0.4 is expressed by color bands. A larger value of $c_{T}$ corresponds to a smaller value of $\tau_{\text{th}}$.

In the same plot, the hadronization time $\tau_{c}$ is also plotted. Combining Eqs. (24) and (52), we obtain

[TABLE]

which is independent of $Q_{s}$ and $\alpha_{s}$. Therefore, $\tau_{c}$ is insensitive to the uncertainty of $c_{T}$ and leads, hence, to the absence of a color band in Fig. 3 for the uncertainty in its value. For the highest RHIC energy, the values of $\tau_{c}$ is only slightly larger than those of $\tau_{\text{th}}$ for $c_{T}=0.18$, with the life time of the thermal QGP phase at most 2 fm/$c$ for the most central collisions. For $c_{T}=0.4$, this is extended to $\simeq 3$ fm/$c$. Clearly, for N_{\text{part}}\raisebox{-2.15277pt}{>\stackrel{{\scriptstyle\textstyle<}}{{\sim}}>}150, the life time of the thermal QGP phase is quite short – the Glasma hadronizes out of equilibrium. For the LHC energy, $\tau_{c}$ is larger by about a factor of two and the life time of the QGP phase is significantly longer.

The temperature at the thermalization time, $T_{\text{th}}=c_{\text{eq}}\,c_{T}\,\alpha_{s}^{2/5}Q_{s}$, is plotted in Fig. 4. As in the case of $\tau_{\text{th}}$, there is a competition between the factors $\alpha_{s}^{2/5}$ and $Q_{s}$. However, the exponent of $\alpha_{s}$ is small here. The variation of $Q_{s}$ wins over that of $\alpha_{s}^{2/5}$ with the result that $T_{\text{th}}$ is an increasing function of $N_{\text{part}}$. In the figure, the crossover temperature $T_{c}=154$ MeV is shown as a black line. For the RHIC energy, the thermal QGP initial temperature for $c_{T}=0.18$ is only slightly larger than $T_{c}$.

IV.3 Comparison of thermal and Glasma photon yields

We now show the comparison of the photon yields in different stages based on the results obtained in Sec. II and III. The thermal photon yield is given by Eq. (25), while the Glasma photon yields in stages (i), (ii) and (iii) are given in Eqs. (36), (45) and (47), respectively. In the expressions for stages (i) and (ii), we set the Coulomb logarithm to be $\mathcal{L}=0.5\log\left(1+2.9/g^{2}\right)$ so that the small-angle approximation reproduces the known formula for the thermal photon yield (see Appendix).

IV.3.1 Photon production rate

To gain insight into photon production in the different stages, we firstly plot, in Fig. 5, the photon production rate $dN/d\tau dy_{p}$ for the most central collisions as a function of time. For simplicity, we fix the value $c_{T}$ to be 0.18. For the Glasma stages (i) and (ii), we plot the the nonequilibrium production rate $\propto\tau I_{g}(\tau)I_{q}(\tau)$ obtained from Eq. (28) (red lines). For the Glasma stage (iii) and the thermal QGP phase, we show the thermal production rate $\propto\tau T^{4}$ (blue solid lines). In the stage (i), the rate decreases as $\tau^{-1}$, while in the stage (ii), the decrease slows down to $\tau^{-1/2}$. Since we do not consider all the possible sources of the photon production, the lines are disconnected between stages (ii) and (iii).131313In Fig. 5, the red and the blue line seem to be connected. However, this is accidental for the current choice of the $Q_{s}^{2}$ profile. As we will discuss below, the Glasma yield strongly depends on the value of $Q_{s}^{2}$, while the thermal yield is insensitive to it. In the stage (iii), the rate increases as $\tau^{5}$ since the temperature increases linearly in time. However, the contribution from the stage (iii) is relatively small because the time duration of this stage is short. In the thermal stage, the rate decreases as $\tau^{-1/3}$. The area under these lines corresponds to the total photon yield $dN/dy_{p}$.

In current hydro simulations for thermal photon production, the bottom-up thermalization scenario is not implemented and hydro modeling is sometimes initialized at early times. For example, in Ref. Paquet:2015lta , the classical Yang-Mills equation is solved with the IP-Glasma initial conditions up to $\tau_{0}=0.4$ fm$/c$ and the system evolution is instantaneously switched to hydrodynamic evolution of a thermal QGP. This situation corresponds to, within our simple model, extending the thermal lines to early times as represented by blue dashed lines in Fig. 5. For comparison, in addition to the bottom-up thermalization scenario, we will consider such hydro scenario extended to the early time, and call the extended hydro stage before $\tau_{\text{th}}$ “early-hydro”.

In the following, we will compare three contributions:

• •

the Glasma contribution in $\tau_{0}<\tau<\tau_{\text{th}}$

• •

the thermal contribution in $\tau_{\text{th}}<\tau<\tau_{c}$

• •

the early-hydro contribution in $\tau_{0}<\tau<\tau_{\text{th}}$.

In the bottom-up thermalization scenario, the total photon yield until the hadronization time is given by the sum of the Glasma contribution and the thermal contribution, while in the hydro scenario that assumes early thermalization, the total yield is the sum of the early-hydro and the thermal contribution.

IV.3.2 Dependence on $Q_{s}$

Thus far, we have fixed the profile of $Q_{s}$, as shown in Fig. 1, by choosing the reference value (the value at the RHIC most central collision) to be 1.4 GeV. If one has a complete description of the space-time evolution in heavy-ion collisions, the effective value of $Q_{s}$ would be fixed for given hadron multiplicity and collision energy. However, as noted previously, because of our lack of knowledge, the reference value of $Q_{s}$ cannot be specified within our framework. We will therefore treat the reference $Q_{s}$ as a free parameter and investigate the dependence of the photon yield on it for given values of the measured charged hadron multiplicity.

We plot in Fig. 6 the bottom-up thermal photon yield (blue solid line) and the Glasma photon yield (red dashed line) as a function of $Q_{s}$. For comparison, the early-hydro photon yield (green dotted line) is also shown in the figure. The values of $N_{\text{part}}$ and $S_{\perp}$ are fixed to those in the most central collisions (centrality 0–5%) and we have used the corresponding experimental data for the charged hadron multiplicity to find the value of the coefficient $c_{\text{eq}}$. For simplicity, we have set $c_{T}=0.18$. We note that the total photon yield until the hadronization time within the bottom-up thermalization scenario corresponds to the sum of the blue and red lines, while the photon yield in the hydro scenario extended to the early time is given by the sum of the blue and green lines. The respective net contributions will be compared later in Fig. 9.

The bottom-up thermal photon yield is not strongly dependent on $Q_{s}$. This can be accounted for by rewriting Eq. (23) as

[TABLE]

where we have used Eq. (52). In this expression, most of the factors are independent of $Q_{s}$. As we discussed previously, $\tau_{c}$ is completely fixed by the measured hadron multiplicity in our model (see Eq. (54)). The coupling $\alpha_{s}$ depends on $Q_{s}$ only logarithmically. Only $\tau_{\text{th}}$ has a power law dependence on $Q_{s}$. However, this dependence is weak due to the competition between the factors $\alpha_{s}^{13/5}$ and $Q_{s}^{-1}$, as noted previously in the context of Fig. 3, thereby explaining the insensitivity of the thermal photon yield to the normalization of $Q_{s}$. The early-hydro photon yield is obtained from (55) by replacing $\tau_{\text{th}}$ by $\tau_{0}$ and $\tau_{c}$ by $\tau_{\text{th}}$. Therefore, it is also insensitive to the value of $Q_{s}$.

In contrast, the Glasma photon yield has a strong dependence on $Q_{s}$. This can be easily understood from the expressions in Eq. (36) and Eq. (45), which are the dominant contributions to the Glasma photon yield. In these expressions for $\frac{1}{Q_{s}^{2}S_{\perp}}\frac{dN}{dy}$, the right hand sides are independent of $Q_{s}$ except for the weak dependence through the coupling. Therefore, $dN/dy$ is approximately proportional to $Q_{s}^{2}$, which is consistent with geometrical scaling of direct-photon production discussed in Klein-Bosing:2014uaa . These results indicate that the preequilibrium Glasma photon production can dominate over the thermal one depending on the value of $Q_{s}$. While the Glasma dominance is pronounced at RHIC, it also begins to dominate at the LHC for $Q_{s}\simeq 3$ GeV.

IV.3.3 Dependence on $N_{\rm part}$

We now turn to the study of the dependence of the photon yield on $N_{\text{part}}$ by fixing the profile for $Q_{s}^{2}$ as shown in Fig. 1. The corresponding values of $Q_{s}$ are indicated in Fig. 6 by black arrows.

In Fig. 7, we plot the thermal photon yield and the Glasma photon yield as a function of $N_{\text{part}}$. Both of the contributions have the uncertainty associated with the undetermined constant $c_{T}$. As for previous figures in Sec. IV.2, this is expressed by color bands corresponding to the range $c_{T}=0.1-0.4$, with the result for the central value $c_{T}=0.18$ denoted by the solid curve.

The RHIC thermal photon yield has a large relative uncertainty. This is because the values of $\tau_{c}$ and $\tau_{\text{th}}$ are very close and the life time of the QGP can be short for a small $c_{T}$ as discussed previously in the context of Fig. 3. On the other hand, the Glasma photon yield has a small uncertainty because the stage (i) and (ii) contributions are independent of $c_{T}$.141414We note that, however, the Glasma contributions involve the other systematic uncertainties, e.g. estimation of the functions $I_{g}$, $I_{q}$ and the Coulomb logarithm $\mathcal{L}$, which are not reflected in the figures. With the exception of the largest values of $c_{T}$, we observe that, in the bottom-up framework, the Glasma contribution is larger than the thermal QGP contribution to the photon yield for the highest RHIC energy. In particular, for off-central collisions, the relative contribution from Glasma becomes more important.

For the LHC energy of $2.76$ TeV, the Glasma photon yield for the most central collisions ranges from $40$–$60$% of the thermal QGP contribution. For $N_{\text{part}}=150$ it is comparable to the central QGP value, and dominates for more peripheral collisions.

As we have seen in Fig. 6, the Glasma photon yield strongly depends on the value of $Q_{s}^{2}$. Therefore the quantitative comparison between the Glasma and the thermal yields also depends on the normalization of the $Q_{s}^{2}$ profile. Nevertheless, for off-central collisions, the qualitative observation that the Glasma contribution is relatively more important is true for any values of $Q_{s}^{2}$.

In Fig. 8, the photon yields in the Glasma stages (i), (ii) and (iii) are plotted separately. For the stage (iii) yield, we fix $c_{T}=0.18$ since the uncertainty is anyway not significant. In stage (i), photons are produced mainly by the scattering among hard gluons and hard quarks, whose transverse momenta are $\sim Q_{s}$. In stage (ii), we have estimated the photon yield by taking into accouunt the scattering between soft gluons and hard quarks. In stage (iii), we evaluated the emission of photons from the bath of soft gluons and quarks. Both of the RHIC and the LHC plots display the same systematics; stages (i) and (ii) give similar contributions to the photon yield and the stage (iii) yield is smaller than the yield from the other two stages. As discussed in Sec. III.3, however, the photon production by the quenching processes of the hard quarks and gluons is not considered in the present study. Such mini-jet photon production processes may have an important contribution in stage (iii).

As we have already noted before, the total photon yield until the hadronization time, within the bottom-up thermalization scenario, is given by the sum of the Glasma and the thermal contributions. On the other hand, the photon yield in the hydro scenario that assumes early thermalization corresponds to the sum of the early-hydro contribution (for $\tau_{0}<\tau<\tau_{\text{th}}$) and the thermal contribution (for $\tau_{\text{th}}<\tau<\tau_{c}$). The photon yields in these two different scenarios are compared in Fig. 9. The contributions that were separately shown in Fig. 6 are now given as sums. The yield in the bottom-up thermalization scenario has the uncertainty associated with the constant $c_{T}$ and it is expressed by color bands. The photon yield in the hydro scenario is naturally independent of $c_{T}$.

The two plots on top of Fig. 9 correspond to the normalization of $Q_{s}$ shown in Fig. 1, which is also indicated in Fig. 6 by black arrows. For this choice of $Q_{s}$ values, at both of the RHIC and LHC energies, the yields in the two scenarios have nearly the same $N_{\text{part}}$ dependence. This agreement is accidental. For larger values of the reference $Q_{s}$, the bottom-up thermalization provides more photons than the hydro model extended to the early time. This is shown in the the bottom plots of Fig. 9, where the normalization of $Q_{s}$ is chosen such that the value at the RHIC most central collision is 2 GeV. The corresponding value at the LHC is 3 GeV.

V Summary and outlook

In this work, we estimated the yields for photon production from both non-equilibrium Glasma stages and the equilibrium thermal QGP stage within the bottom-up thermalization scenario of heavy-ion collisions. While the uncertainties from our lack of knowledge of the coefficients multiplying parametric estimates are large, they can be constrained significantly by the measured charged hadron multiplicity. For an assumed $Q_{s}^{2}$ profile as a function of centrality (Fig. 1), we found that at RHIC energies the Glasma photon yields are larger than the thermal photon yields for a wide parameter range. This dominance is especially pronounced for more peripheral collisions though we must caution that weak coupling estimates for RHIC energies are likely not reliable for peripheral collisions. At the LHC, the thermal photon yields are larger, but even at the most central collisions the Glasma contribution can range from 40–60% of the thermal QGP yield. If the reference $Q_{s}$ in the $Q_{s}^{2}$ profile is increased, we find that the Glasma contribution become larger than the QGP photon yield even for central collisions.

We also made a comparison between the photon yields in the bottom-up thermalization scenario to those in a hydro scenario that assumes the system thermalizes at the early time $\simeq Q_{s}^{-1}$. For the $Q_{s}^{2}$ profile shown in Fig. 1, the two scenarios give comparable photon yields. If we assume a larger value (by about 50%) of the reference $Q_{s}$, the bottom-up thermalization scenario provides a larger photon yield relative to the early-hydro scenario.

Our results point to the urgent need for more refined computations to reduce the uncertainties we identified in the computation of photon yields from different stages of the Glasma evolution in addition to more sophisticated computations of thermal photon yields. These include first principles classical-statistical computations of the photon yields in the first stage of Glasma evolution as well as more detailed kinetic theory computations that match to these and to viscous hydrodynamics at later times.

We have only considered photon yields in this work. Computations of photon spectra and flow coefficients within the bottom-up framework, and their comparison to the available data, can help either rule out the bottom-up framework of QGP equilibration or at least strongly constrain the viable set of free parameters. The possibility that these can be extracted with increasing precision offers the promise that a quantitative theory of the equilibration of strongly correlated quark and gluon matter can be developed further and tested.

Acknowledgements

We would like to thank Oscar Garcia-Montero, Niklas Mueller, Chun Shen, Bjoern Schenke, Soeren Schlichting and Prithwish Tribedy for very valuable discussions and comments. R. V. is supported under DOE Contract No. DE-SC0012704. He would like to thank the Institut für Theoretische Physik, Universität Heidelberg for kind hospitality and support via the Excellence Initiative during the early stages of this work. This work is part of and supported by the DFG Collaborative Research Centre “SFB 1225 (ISOQUANT)”.

Appendix A Small-angle approximation for the photon production

The scattering amplitude due to the exchange of a massless particle has an infrared divergence when the exchanged momentum goes to zero. Indeed, the amplitude for the pair annihilation process in Eq. (3) and that of Compton scattering in Eq. (4) diverge for $t\to 0$ or $u\to 0$. The small-angle approximation is applicable to such collision processes landau1936kinetic ; lifshitz1981physical . When an incoming particle and an outgoing particle that is kicked by small momentum exchange are of the same species, the collision integral for that process can be approximated as a diffusion term. When the incoming and the outgoing particle are different species, the collision integral is replaced by a source term, which has a simple form Blaizot:2014jna . Photon production by the pair annihilation and Compton scattering corresponds to this latter case.

Following the same procedure as the calculation outlined in Blaizot:2014jna , we apply the small-angle approximation to the photon production formula in Eq. (1). First, let us consider the pair annihilation process. Since the $t$-channel and the $u$-channel give the same contribution, it is sufficient to compute the $t$-channel contribution alone, multiplied by a factor of two. This gives

[TABLE]

We expand kinematic variables in terms of the exchanged momentum ${\boldsymbol{q}}={\boldsymbol{p}}-{\boldsymbol{p}}_{1}$ (henceforth asssumed to be small in magnitude) to obtain,

[TABLE]

with ${\boldsymbol{v}}_{1,2}={\boldsymbol{p}}_{1,2}/p_{1,2}$. Furthermore,

[TABLE]

Keeping the leading order terms in $q$, one obtains

[TABLE]

where, similarly to ${\boldsymbol{v}}_{1}$ and ${\boldsymbol{v}}_{2}$, ${\boldsymbol{v}}={\boldsymbol{p}}/p$. By a straightforward computation, one can show that the ${\boldsymbol{q}}$-integration is independent of ${\boldsymbol{v}}$ and ${\boldsymbol{v}}_{2}$, and express it as

[TABLE]

where the logarithmic divergence can be further expressed as

[TABLE]

In thermal field theory, the IR cutoff is given by the Debye mass scale $m_{D}\sim g^{2}T^{2}$ and the ultraviolet cutoff is given by the temperature $T$. Hence $\mathcal{L}\sim\log(1/g)$.

The production rate from the annihilation process simplifies in the small-angle approximation to

[TABLE]

The expression within the integrand corresponds to the density of scatterers (quarks/anti-quarks) enhanced by the Bose factor of the final state gluons.

For the Compton scattering contribution, we can neglect the $s$-channel contribution in this approximation. Following the same procedure as in that for the annihilation process, one can derive,

[TABLE]

The expression within the integrand in this case corresponds to the density of scatterers (gluons) suppressed by the Pauli factor of the final state quarks/anti-quarks. Summing the Compton and annihilation contributions, we obtain

[TABLE]

Since the small-angle approximation computation of the photon rate only involves kinematic approximations based on the dominant contributions to on-shell $2\leftrightarrow 2$ quark-gluon scattering, our result is applicable to either equilibrium or non-equilibrium situations where kinetic theory is applicable. To check the validity of this approximation, let us consider thermal equilibrium and compare our approximate result in Eq. (68) with the thermal photon rate we quoted in Eq. (6). For the equilibrium distribution (with vanishing chemical potential),

[TABLE]

The thermal rate with the small-angle approximation then gives

[TABLE]

If we identify the Coulomb logarithm $\mathcal{L}$ with the logarithmic factor in Eq. (6), two results agree up to a numerical factor of 2. We note, however, that this discrepancy can be traded for the uncertainty in ${\mathcal{L}}$.