A Lattice Calculation of Thermal Dilepton Rates

Using clover improved Wilson fermions we calculate thermal vector meson correlation functions above the deconfinement phase transition of quenched QCD. At temperatures 1.5 Tc and 3Tc they are found to differ by less than 15% from that of a freely propagating quark anti-quark pair. This puts severe constraints on the dilepton production rate and in particular rules out a strong divergence of the dilepton rate at low energies. The vector spectral function, which has been reconstructed using the Maximum Entropy Method, yields an enhancement of the dilepton rate over the Born rate of at most a factor two in the energy interval 4<E/T<8 and suggests that the spectrum is cut-off at low energies by a thermal mass threshold of about (2-3)T.


Introduction
The dilepton spectrum is one of the key observables in the study of thermal properties of the hot and dense medium created in heavy ion collisions [1,2,3]. Thermal modifications of the dilepton spectrum have been observed at large invariant masses where they reflect the suppression of heavy quark resonances [4]. At low energies it is expected that non-perturbative in-medium modifications of the quark anti-quark interactions as well as quark anti-quark annihilation in the quark gluon plasma influence the thermal dilepton rate. Perturbative 2-loop [5] calculations as well as calculations based on hard thermal loop (HTL) resummed perturbation theory [6] lead to an enhancement of the spectrum at low energies, which dominates over the suppression arising from massive quasi-particle poles [1,6] or non-perturbative gluon condensates [7]. Indeed, an enhancement of low mass dilepton production has been observed in heavy ion experiments [8].
The thermal dilepton rate describing the production of lepton pairs with energy ω and total momentum p is related to the Euclidean correlation function of the vector current, J µ V (τ, x) =ψ(τ, x)γ µ ψ(τ, x), which can be calculated numerically in the framework of lattice regularized QCD. This relation is established through the vector spectral function, σ V (ω, p, T ), which directly gives the differential dilepton rate in two-flavor QCD, and at the same time is related to the Euclidean vector meson correlation function, where the Euclidean time τ is restricted to the interval [0, 1/T ]. A direct calculation of the differential dilepton rate thus becomes possible, if the above integral equation can be inverted to determine σ V (ω, p, T ). Although finite temperature lattice calculations, which are performed on lattices with finite temporal extent N τ , will generally provide information on G V (τ, p, T ) only for a discrete and finite set of Euclidean times τ T = k/N τ , k = 1, ..., N τ , this may be achieved through an application of the maximum entropy method (MEM) [9,10]. We will present first results on such a calculation in this letter. However, even without invoking statistical tools like MEM accurate numerical results on G V (τ, p, T ) will themselves provide stringent constraints on spectral functions calculated in other perturbative  Table 1: Simulation parameters for the calculation of G V (τ, T ) on lattices of size 48 3 × 12 and 64 3 × 16. The third column gives the gauge coupling β = 6/g 2 . The values used for the clover coefficient c sw , the critical hopping parameter κ c and the current renormalization constant Z V are discussed in the text. The values for T /T c given in the first column are based on results for the non-perturbative β-function given in Ref. [14].
or non-perturbative approaches. This in turn will put constraints on the thermal dilepton rates.
We present here results of a calculation of vector current correlation functions and reconstructed spectral functions at two temperatures above T c , i.e. T = 1.5 T c and 3 T c . We will restrict our discussion to the zero momentum limit ( p = 0) of these quantities and we will therefore suppress any explicit momentum dependence in the arguments of G V and σ V in the following.

Thermal correlation functions
The correlation functions G V (τ, T ) have been calculated within the quenched approximation of QCD using non-perturbatively improved clover fermions [12,13]. This eliminates the O(a) discretization errors in the fermion sector. Moreover, in order to explicitly control the cut-off dependence of the numerical results we have used lattices of size N 3 σ × N τ with fixed ratio of spatial to temporal lattice extent, N σ /N τ = 4, and two different temporal lattice sizes, N τ = 12 and 16, respectively. Due to the large temporal extent of the lattice, the lattice spacing is rather small already at the deconfinement transition and becomes even smaller for larger temperatures. Using T c = 265 MeV [14] to set a scale for the lattice spacing one finds, a = 1/N τ T ≃ (0.75T c /T N τ ) fm. For the two different temperatures and lattice sizes used here it ranges from a = 0.015 fm to a = 0.04 fm. The comparison of results obtained on lattices with different temporal extent and fixed T indeed confirms that the results for G V (τ, T ) are not significantly influenced by finite cut-off effects.
The thermal correlation functions have been constructed from quark propagators obtained from the inversion of the clover improved Wilson fermion matrix. As we τT  have performed calculations in the high temperature, chirally symmetric phase of QCD there are no massless Goldstone modes which could interfere with the convergence of the iterative solvers for this matrix inversion. We thus could perform calculations directly in the limit of vanishing quark masses, i.e. directly at the critical values of the hopping parameter κ c given in Tab. 1. The values for c sw as well as critical κ values were obtained from [13] and interpolation. Note that a finite temperature critical κ defined by a vanishing quark mass might differ from the T = 0 value κ c by finite lattice spacing and finite volume effects. Indeed, the values for the quark mass obtained from the PCAC relation [15] are very small but not exactly zero, see Tab. 1. However, since the correlation functions above T c are almost quark mass independent close to the chiral limit [16] it is not crucial to hit κ c (T = 0) precisely.
Some care has to be taken over the proper renormalization of the vector current J V µ used in the lattice calculations. In order to be able to compare results obtained with different N τ and T as well as with calculations performed in different continuum regularization schemes we use the renormalized local current J renorm Here Z V (g 2 ) is the current renormalization constant for which we use the nonperturbative values determined in [13]. The overall error on the normalization [13,17], including the effect of the quark mass not being exactly zero has been estimated to be less than 1 %.
In the high temperature limit the ratio G V /T 3 can be calculated perturbatively. The leading order result is obtained by using free, massless quark propagators in the calculation of the current-current correlation function. For vanishing momentum of the quark anti-quark pair ( p = 0) the corresponding spectral function is given by σ V (ω, T ) = 0.75 π −2 ω 2 tanh(ω/4T ). In this limit the differential dilepton rate reduces to the Born rate, In Fig. 1a we show our numerical results for G V /T 3 calculated at T /T c = 1.5 and 3. These results are compared to the lattice version of the free vector correlation function, G free,lat V /T 3 [18]. It is evident that the correlators obtained at two different temperatures agree with each other quite well and stay close to the leading order perturbative (free) vector correlation function already quite close to T c . The deviations from G free,lat V /T 3 are better visible in the ratio G V /G free,lat V which is shown in Fig. 1b. It is noteworthy that G V (τ )/T 3 is definitely larger than the corresponding free correlator, G free,lat V /T 3 for both temperatures and all Euclidean times 0 < τ T < 1. This is furthermore supported by a calculation of these ratios for different values of the cut-off, e.g. for N τ = 12 and 16 at the same value of T /T c . Results agree within statistical errors.
As is obvious from Eq. 2 the correlation function G V (τ, T ) receives contributions from all energies ω. Due to the presence of the integration kernel K(τ, ω) = ch(ω(τ − 1/2T ))/sh(ω/2T ) large energies are, however, exponentially suppressed. This suppression becomes more efficient with increasing Euclidean time τ T . For τ T ≃ 1/2 the correlation function thus is most sensitive to the low energy part of σ V (ω, T ). In fact, using the free spectral function in Eq. 2 we estimate that energies larger than ω/T ≃ 15 contribute only 2% to the central value G V (1/2T, T )/T 3 whereas in the low energy regime, ω/T < ∼ 3, the corresponding contribution still amounts to 12%. G V (1/2T, T )/T 3 thus is sensitive to modifications of the vector spectral functions up to ω/T ≃ 15 and provides a stringent test of models and approximate calculations of the low energy part of the vector spectral function. We find which is about 10% larger than the free, infinite volume value G free V (1/2T, T ) ≡ 2 as well as the corresponding free value on a finite 64 3 × 16 lattice, G free,lat V (1/2T, T ) ≡ 2.0368. Our result thus suggests that at least for a significant range of energies the vector spectral function is enhanced over the free case. Such an increase cannot easily be incorporated in simple quasi-particle pictures, which only include massive quasiparticle poles in the quark spectral functions. They lead to a reduction of G V (τ, T ) relative to G free V (τ, T ) and in turn to a reduced differential dilepton rate [1,7,19]. On the other hand, HTL-resummed perturbative calculations, which incorporate not only poles at non-zero quasi-particle masses but also additional contributions from cuts in the HTL-resummed quark spectral functions, show that both nonperturbative features of quarks propagating in a hot plasma play an important role for the low energy part of the meson spectral functions [6]. The cut contributions in the HTL approach [6] as well as in 2-loop perturbative calculations [5] even lead to divergent vector spectral functions at small energies which also makes G V (τ, T ) to diverge for all Euclidean times τ . Assuming that at our aspect ratio of T V 1/3 = 4 finite volume effects do not play an important role, the mere fact that we do find a finite result for the vector correlator via Eq. 2 implies that σ V (ω, T ) vanishes in the limit ω → 0.

Reconstructed spectral functions
The correlation functions G V shown in Fig. 1 clearly stay close to the leading order perturbative result. This suggests that also the spectral functions σ V are close to that of the free case. We have reconstructed σ V from G V using MEM, which has successfully been applied at T = 0 [9,11] and also has been tested at finite temperature [10]. In order to take into account finite lattice cut-off effects in the reconstruction of the spectral function we introduce a lattice version, K L (τ, ω, N τ ), of the continuum kernel K(τ, ω) which appears in Eq. 2. With the lattice kernel the spectral function σ V for vanishing momentum is defined through where K L (τ, ω, N τ ) is the finite lattice approximation of a kernel appropriate to describe the correlation function of a free boson at finite temperature, This lattice kernel differs from the T = 0 version used in [9] in so far as we have explicitly taken into account the finite temporal extent of the lattice. Of course, K L approaches K in the limit of large N τ . Moreover, for the energy range important for the description of the correlation functions close to τ T = 1/2, i.e. for ω/T < 15, they differ by less than 2% already on lattices with temporal extent N τ = 16. Using this lattice kernel we have checked that the MEM analysis of correlation functions constructed from given input spectral functions successfully reproduces these input functions. In particular, we find that the free vector spectral function for massless In Fig. 2a we show spectral functions reconstructed from the vector current correlation functions at T /T c = 1.5 and 3 on the 64 3 × 16 lattices. In the MEM analysis we have taken into account energies up to ω/T = 4N τ and we used the lowest order zero temperature perturbative result m(ω) = 0.75 π −2 ω 2 as a default model. In order to test the statistical significance of our results we have split our total data sets in 8 jackknife blocks and performed a MEM analysis on each of these blocks. The resulting jackknife error, which is given as a band in the insertion of Fig. 2a, shows that our data sample is large enough to yield statistically significant results. We also have checked that a variation of the default model by 20% only leads to minor changes in the low energy part of σ V , which stay within the error band shown in Fig. 2a. Another question is to what extent the spectral function obtained from the MEM analysis is unique. As the MEM analysis is based on a finite data set the reconstruction of the spectral function itself is not unique. This uncertainty is incorporated in the MEM-error, which is calculated from the covariance matrix of the spectral function [20] for four energy bins in the interval 0 ≤ ω/T ≤ 16. The resulting error on the average value of σ V (ω, T )/ω 2 in these bins is also shown in Fig. 2a.
The reconstructed spectral functions show a broad enhancement over the free spectral function for ω/T > ∼ 16, i.e. ωa > ∼ 1. This regime clearly is influenced by lattice cut-off effects. In fact, a similar enhancement is observed in T = 0 spectral functions and has been attributed to unphysical contributions of the heavy fermion doublers present in the Wilson (and clover) lattice fermion action [11]. As discussed above this regime of large energies does not contribute to the central values of the vector correlator for τ T ≃ 1/2. The relevant energy regime is given in the insertion of Fig. 2a. The enhancement observed for G V (1/2T, T ) thus is due to the peak in σ V (ω, T )/ω 2 found for ω/T ≃ (5 − 6). We furthermore note that for both temperatures the spectral functions have a sharp cut-off at small energies. They drop rapidly below ω/T ≃ 5 and vanish below ω/T ≃ 3. This is in contrast to perturbative calculations of σ V (ω, T ), which lead to a divergent spectral function for small energies.
In Fig. 2b we show the differential dilepton rate calculated from Eq. 1 and the low energy part of σ V shown in Fig. 2a. The comparison with the Born rate shows that for all energies ω/T > ∼ 4 the difference is less than a factor 2. For energies ω/T < ∼ 3 the dilepton rate drops rapidly and reflects the sharp cut-off found in the reconstructed spectral functions. Of course, with our present analysis, which is based on rather small statistics, we cannot rule out the existence of sharp resonances like the van Hove singularities present in the HTL-resummed calculations. It also may well be that the broad peak found by us will sharpen with increasing statistics as it has been observed in related MEM studies of hadron correlation functions at zero temperature [9]. Nonetheless, any further sharpening of resonance like peaks in a certain energy interval of the spectral function has to be compensated by an even closer agreement with the free spectral function in other energy intervals due to the constraint given by Eq. 4.

Conclusions
To conclude, we find that already close to the QCD phase transition temperature, i.e. for T = 1.5T c and 3T c , the vector correlation function deviates only by less than 15% from the leading order perturbative result. This also is reflected in the reconstructed spectral function which deviates by less than a factor two from the leading order perturbative form for energies ω > ∼ 4T . The most pronounced feature of the spectral function and the resulting dilepton rate is the presence of a sharp cut-off at low energies, ω ∼ (2 − 3)T . If a threshold of similar magnitude persist also closer to T c , there will be no thermal contribution to the dilepton rate at energies ω < ∼ 2T c during the entire expansion of the hot medium created in a heavy ion collision. This is consistent with present findings at SPS energies [8] and may become visible as a threshold in the dilepton rates for long-lived plasma states expected to be created at RHIC or LHC energies.