Study of Quarkonium in QGP from Unquenched Lattice QCD

This paper discusses the charmonium and bottomonium correlators in the pseudoscalar channel and the corresponding spectral reconstruction on the lattice. The absence of a transport peak in the pseudoscalar channel spectral function allows for an easier study of the in-medium modification of bound states. However, extracting spectral information from Euclidean correlators is still a numerically ill-posed problem. To constrain the spectral reconstruction, we use an ansatz motivated from perturbation theory. The perturbative model spectral function has two main contributions: a thermal part around the threshold obtained from pNRQCD and the vacuum part well above the threshold. These two regions are matched continuously, and the model spectral function is obtained by introducing parameters that control the overall thermal shift of the peak and the overall amplitude. The lattice correlator data is computed using clover-improved Wilson valence fermions on large and fine gauge field configurations generated using Nf=2+1\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$N_f=2+1$$\end{document} flavors Highly Improved Staggered Quark action with physical strange quark mass ms\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$m_s$$\end{document}, and slightly heavy degenerate up and down quark masses ml=ms/5\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$m_l=m_s/5$$\end{document} that correspond to mπ≃320\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$m_\pi \simeq 320$$\end{document} MeV. Our results obtained at T=220\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$T=220$$\end{document} MeV and T=251\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$T=251$$\end{document} MeV suggest that no resonance peaks are needed to describe the charmonium lattice data at these temperatures, while for bottomonium thermally broadened resonance peaks persists.


Introduction
Heavy quark-antiquark bound states, quarkonia, like charmonium and bottomonium, serve as a good thermometer for the quark-gluon-plasma (QGP) created in ultrarelativistic heavy ion (AA) collisions [1].This is due to the fact that they are produced during the early stages of the collisions and participate the entire evolution of the QGP.Some of them will dissociate into heavy quarks when experiencing temperature increasing while some can remain as bound states due to their hierarchically small sizes and large binding energies [2,3].How to theoretically understand the suppression of quarkonia yields observed in AA collisions compared to that in the proton-proton (pp) collisions [4][5][6], however, remains unclear.This is mainly due to the complicated entanglement between the cold and hot nuclear effects [7].Given that the suppression resides in the non-perturbative regime, naturally we want to interpret it using lattice QCD calculations.Besides, the experiments at RHIC and LHC reveal that open heavy mesons show an unexpectedly substantial elliptic flow, comparable to that of light mesons [8][9][10].This indicates that heavy quarks flow as efficiently as light quarks do [11][12][13].Trying to interpret this requires a modeling of the heavy quark diffusion in a hot and dense medium, which is challenging in perturbation theory [14][15][16][17], but feasible from lattice QCD.Studies carried out in the heavy quark mass limit can be found in [18][19][20][21][22][23] in both quenched approximation and (2+1)-flavor lattice QCD.Accessing the first mass suppressed correction to the diffusion coefficient are also attempted [22,[24][25][26].
Tackling the above puzzles relies on the knowledge about the spectral function of quarkonia, which can be extracted from the fully relativistic Euclidean mesonic correlation functions.The spectral function in the vector channel contains all information about the in-medium hadron properties like the dissociation temperatures and heavy quark diffusion, and was intensely studied in the past decade, see [27][28][29][30][31][32][33] for a selection.While the one in the pseudoscalar channel provides an easier access to the fate of the bound states due to the vanishing of the transport peak [34,35].This was investigated in the quenched case in [36].A review of current status of lattice studies on heavy quarks and quarkonium in extreme conditions can be found in [37][38][39].
In this work we aim to provide a (2+1)-flavor lattice QCD calculation for the charmonium and bottomonium correlation functions in pseudoscalar channel, G PS (τ ), at two temperatures T =220 MeV and T =251 MeV.The Euclidean correlator is related Ref. [40].The tuning of the quark masses is done at lowest temperature with maximum temporal extent Nτ =64.
to the spectral function via following equation Due to statistical errors in G PS (τ ) and a finite resolution of the lattice, the solution of ρ PS is not unique.However, the spectral reconstruction can be constrained with further inputs.In this study we model the spectral function with a perturbatively inspired ansatz, interpolating between the thermal part around the threshold obtained from pNRQCD [41] and the vacuum part well above the threshold [36].This article is organized as follows.We begin with the lattice setup and clover mass tuning, which are outlined in Sec.(2).Sec.(3) presents the numerical data of lattice correlators, while Sec.(4) details the process of obtaining perturbative spectral functions.In addition, Sec. ( 5) provides an explanation of the spectral reconstruction procedure.Finally, we summarize our findings and present future outlook in the concluding section.

Lattice parameters and mass tuning
Quarkonium correlators are obtained at various temperatures, as specified in Tab. 1, using clover-improved Wilson fermions as valence quarks.The lattice correlation functions are computed on large and fine gauge field configurations generated with N f = 2 + 1 and with the Highly Improved Staggered Quark (HISQ) action [42] and tree-level improved Lüscher-Weisz gauge action [43,44].The sea quark masses are tuned to the physical strange quark mass, m s , and degenerate up, down quark masses m l = m s /5, corresponding to m π 320 MeV.The fermion valence action is tadpoleimproved by incorporating clover-improvement, which reduces the O(a) cut-off effects.To this aim, the Sheikholeslami-Wohlert coefficient [45] is chosen as c SW = 1/u 3 0 , where u 0 is the tadpole factor determined by taking the fourth root of the plaquette expectation value calculated on our HISQ lattices.
The masses of the lattice calculated correlators are determined by the hopping parameter κ, which should be tuned to match the quarkonium meson mass spectrum observed in experiments [46].In the heavy quark sector, the tuning is performed with respect to the spin-averaged charmonium m cc = (m ηc + 3m J/ψ )/4 and bottomonium represents the physical mass in lattice units.The values of κ that are tuned for charm and bottom are 0.13164 and 0.11684, respectively.

Correlators
We present the pseudoscalar channel charmonium and bottomonium correlators in Fig. 2. Instead of the correlator itself, we used the ratio G PS (τ )/G free PS (τ ) for a better visibility of the behavior as a function of τ at various temperatures, where G free PS (τ ) reads [36] where m(μ ref ) denotes the running mass at reference scale μref =2 GeV.We choose the values of M q to be 1.5 GeV for charmonium and 4.7 GeV for bottomonium.One can see that the charmonium suffers more thermal effects than the bottomonium does, because the charmonium correlators at high temperatures deviate more from the low-temperature ones.

Perturbative spectral functions in PS-channel
To construct the quarkonium spectral functions in the pseudoscalar channel we follow the strategy developed in the quenched approximation [36].For ω 2M , the spectral function is insensitive to temperature effects, therefore it is replaced by its vacuum counterpart.The vacuum spectral function, normalised by ω 2 m 2 (μ), in the Fig. 2 The ratio of pseudoscalar lattice correlator with its corresponding free correlator at different temperatures for charmonium (left) and bottomonium (right).Statistical errors are estimated using the Jackknife procedure.
pseudoscalar channel reads Here R-function is valid up to O(α 3 s ) and it is defined in Ref. [36].One can also express the spectral function ρ PS normalized by the pole mass M as follows: where R is again given in Ref. [36].However, this expansion does not converge as the correction term becomes increasingly large as ω increases.In contrast, the form in Eq. ( 3), obtained by re-expanding the pole mass in terms of the M S mass in Eq. ( 4), leads to smaller corrections at large ω .Following Ref. [36] we can then estimate the pole mass self-consistently by equating Eq. (3) and Eq. ( 4), which is given by, Here the running of the M S mass has been calculated using 5-loop running, starting from m(μ = 2 GeV).Tab. 2 shows the values of M x obtained at O(α 2 s ).Now we consider thermal part of the spectral function around the threshold i.e. ω ≈ 2M .In the heavy quark limit, relativistic effects in quark-antiquark bound states become negligible, allowing for a non-relativistic potential description.For separations much smaller than 1/gT , the potential is effectively a zero temperature potential.However, as the separation approaches ∼ 1/gT , the contribution of soft gluon exchange  2 The quark masses relevant for the perturbative spectral functions for full QCD, utilizing 5-loop running for α s and setting Λ MS to 0.339 GeV [50].The mass Mx was defined using Eq. ( 5), with x = 4...8.For the charmonium case, we selected m(μ ref ) = 1.08 GeV to ensure that m(μ = m) yields mass in the MS scheme that is relevant.
On the other hand, a value of m(μ ref GeV is close to the value relevant for the bottomonium case.
becomes significant, requiring a Hard Thermal Loop resummation.As a result, the potential in this region becomes temperature dependent and acquires a complex component.In this non-relativistic regime, one can get the pseudoscalar spectral function from the vector channel as The vector channel spectral function is given by where C > is a Wightman function and it is obtained by solving the following Schrödinger equation with initial conditions C V > (0; r, r ) = 6N c δ (3) ( r − r ) .The potential V T (r) for positive t has the following form [47][48][49] Here φ(x) is defined to be Here α S is the thermal coupling and m D (∼ gT ) is the Debye mass.Physically, the real part of the potential is related to Debye screening in the plasma.It affects the thermal mass shift of bound states.On the other hand, the imaginary part originates from Landau damping of space like gauge fields, which causes the bound state peak to broaden.At a much lower ω than the threshold, the Schrödinger description described above overestimates the spectral function, as this formalism is no longer valid in this region.
A naive NLO calculation shown in Ref. [36] predicts that the spectral function in this region is exponentially suppressed.To model this suppression within the Schrödinger formalism itself, we multiply the imaginary part of the potential by e −|ω−2 M |/T for ω < 2 M .Moreover, at small spatial distances, r 1/m D , the thermal potential is replaced by a vacuum one [36].Fig. 3 The comparison between perturbative and model spectral functions.The dotted lines represent the perturbative spectral functions, while the solid lines illustrate their modifications as per Eq. ( 12), for both charmonium (on the left) and bottomonium (on the right).The solid curves correspond to the imaginary-time correlators, labeled as "x=mod" in Fig. 4, which exhibit good agreement with the lattice data for all distances, except at small distances where the lattice cut-off effects dominate due to the finite lattice spacing.
After computing both thermal and vacuum contributions to the spectral function, we combine them as In Eq. ( 11) we introduce a multiplicative factor A match to the thermal part and it is determined such that both thermal and vacuum parts of the spectral function are connected smoothly at some ω = ω match , where 2M < ω match < 3M .The resulting pseudoscalar spectral functions at various temperatures for both charmonium and bottomonium are shown in Fig. 3. Unlike charmonium, the bottomonium has a resonance peak that is visible at temperatures 220 MeV and 251 MeV.The resonance peak for bottomonium hints the existence of bound state at these given temperatures, however the peak might get broaden at higher temperatures.3 Estimates of the best fit parameters at two temperatures according to Eq. ( 12) for charmonium and bottomonium.

Modelling and spectral reconstruction
In this section, we discuss the process of modeling perturbative spectral functions.We introduce two parameters, A and B, in the perturbative spectral function shown in Fig. 3. Parameter A will account for the normalization of the correlator data, while parameter B will allow for the adjustment of thermal mass shift because pole masses are poorly determined in perturbation theory.As a results the model spectral function gets the following form: Before fitting the lattice data with the above ansatz, it is important to note that the pseudoscalar correlator can be multiplicatively renormalized, which only affects the fit parameter A. Since our analysis involves correlators at different temperatures obtained from a single lattice spacing, the renormalization constant will be the same for all temperatures.Therefore, we do not need to renormalize the lattice correlator for our present study.Nevertheless, to avoid obtaining large values of A after fitting, we normalize the lattice correlator such that it matches the perturbative correlator at a very short distance of τ T = 0.107 at a temperature T = 220 MeV.We perform uncorrelated χ 2 fit of this normalized lattice data with the ansatz in Eq. (12).The short distance part of the lattice correlator are largely affected by lattice artifacts.As a result, we fit the data starting from τ T =0.25 for T = 220 MeV and τ T =0.21 for T = 251 MeV.The resulting fit parameters, A and B, for two distinct temperatures, are presented in Tab. 3. Preliminary results of the model spectral function ρ mod and the lattice correlator data, are shown in Fig. 4. The important observation is the shift of the thresholds towards larger energy in the model spectral function.Overall, the observed shift is perfectly admissible and could be due to the uncertainties in the determination of the pole mass or the thermal mass corrections.It is important to take these uncertainties into account when interpreting the results.

Conclusion and Outlook
In this work we have presented preliminary results of the spectral reconstruction in the pseudoscalar channel for quarkonium, using an ansatz for the fitting procedure based on perturbation theory.We solved the Schrödinger equation for the thermal part of the spectral function around the threshold, using the perturbative finite temperature potential.For energies well above the threshold, we considered the vacuum part and matched both parts smoothly.Our analysis was limited to two temperatures, namely T = 220 MeV and T = 251 MeV, and we obtained lattice correlators using clover-improved Wilson fermions measured on large N f = 2 + 1 HISQ gauge field configurations generated at gauge coupling β = 8.249.we have tuned the quark masses so that the experimental meson spectrum was reproduced on the lattice within errors.
Fitting the spectral function model to the lattice data, we observe a good description of the lattice correlators suggesting that no resonance peaks are needed to describe the charmonium lattice data at these temperatures, while for bottomonium thermally broadened resonance peaks persist at the analyzed temperatures.While our model describes the lattice data well for large time-slice separations, we observed a discrepancy at small distances, which probably is due to the cut-off effects.To estimate these effects, we are in the process of adding more lattice spacings followed by a continuum extrapolation.The non-perturbative effects are important in the temperature range we are considering in this work.Therefore, to improve the model of the spectral function for this range in the future, we will estimate the thermal part of the spectral function using non-perturbatively determined thermal potential from the lattice [52][53][54].Additionally, we plan to estimate the pole mass, which affects the thermal part of the spectral function and is responsible for the peak position, from the zero-temperature Cornell potential [51].This will help us better understand the thermal mass shift.
Furthermore, we plan to include lattices at physical quark masses and to consider the vector channel of the quarkonium spectral function, from which we will be able to estimate the heavy quark diffusion coefficients.
to Luis Altenkort for generating gauge field configurations and integrating meson measurements into the QUDA code.This work is supported by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation)-Project number 315477589-TRR 211.The computations in this work were performed on the GPU cluster at Bielefeld University using SIMULATeQCD suite [55,56] and QUDA [57].We thank the Bielefeld HPC.NRW team for their support.

Fig. 1 Fig. 1
Fig.1Charm quark mass tuning on the mixed action.The solid line represents the experimental value of the spin-averaged charmonium in lattice units taken from the PDG[46].Whereas triangles represent masses from lattice calculations for different values of κ in lattice units.

Fig. 4
Fig.4The data points correspond to the lattice correlators, while the lines depict the model correlators obtained from Eq. (1) for charmonium (on the left) and bottomonium (on the right).

Table 1
Lattice parameters for N f =2+1, m l = ms/5 HISQ configurations for different temperatures.The lattice spacing a at β=8.249 is obtained using f K -scale parametrization, cf.