Static quark-antiquark potential in the quark-gluon plasma from lattice QCD

Static quark-antiquark potential in the quark-gluon plasma from lattice QCD

Yannis Burnier Institute of Theoretical Physics, EPFL, CH-1015 Lausanne, Switzerland    Olaf Kaczmarek Fakultät für Physik, Universität Bielefeld, D-33615 Bielefeld, Germany    Alexander Rothkopf Institute of Theoretical Physics, Universität Heidelberg, Philosophenweg 12, D-69120 Germany
July 14, 2019

We present a state-of-the-art determination of the complex valued static quark-antiquark potential at phenomenologically relevant temperatures around the deconfinement phase transition. Its values are obtained from non-perturbative lattice QCD simulations using spectral functions extracted via a novel Bayesian inference prescription. We find that the real part, both in a gluonic medium as well as in realistic QCD with light , and quarks, lies close to the color singlet free energies in Coulomb gauge and shows Debye screening above the (pseudo) critical temperature . The imaginary part is estimated in the gluonic medium, where we find that it is of the same order of magnitude as in hard-thermal loop resummed perturbation theory in the deconfined phase.


The potential acting between a heavy quark and anti-quark in a thermal medium is a central ingredient in our understanding of the strong interactions, described by quantum chromo-dynamics (QCD). The bound states it sustains, heavy quarkonium, are precision probes connecting theory and experiment Brambilla:2010cs (). They allow us to test QCD via low temperature spectroscopy Gray:2005ur (), as well as through their in-medium modification Shuryak:1978ij (); Asakawa:2003re (); Jakovac:2006sf () observed in the quark gluon plasma created in relativistic heavy ion collisions. In particular the open question of melting and regeneration observed at RHIC and LHC Adare:2006ns () urges a quantitative understanding of their in-medium behavior.

A wealth of intuition has been accumulated in the past based, in part, on analogies with Abelian theories Shuryak:1978ij (), potential modeling Nadkarni:1986as () and strong coupling approaches SC (). Lattice QCD at tells us Koma:2006si () that the potential rises linearly before flattening off due to string breaking. Perturbation theory on the other hand shows that Debye screening plays a major role in the deconfined phase. At , reached in current experiments, we expect that the medium gradually weakens the interaction. How the transition between the two regimes manifests itself quantitatively in the potential however remained unanswered. Due to recent conceptual and methods developments we are now able to present in this letter a first principles determination of the temperature dependence of the static inter-quark potential in the phenomenologically relevant, i.e. non-perturbative regime around the phase transition.

The advent Barchielli:1986zs () of modern effective field theory allowed to put the definition of the static potential on a rigorous mathematical footing. By exploiting the separation between the heavy quark rest mass and medium scales, a derivation from a dynamical QCD observable, the real-time thermal Wilson loop was achieved,


This expression has been evaluated at finite temperature in hard thermal loop (HTL) resummed perturbation theory Laine:2007qy () and was found to be complex valued. In the deconfined phase the real part shows Debye screening, while the imaginary part is related to the scattering (Landau damping) and absorption (singlet-octet transition) of gluons from the medium. Even though at leading order the real part coincides with the color singlet free energies in Coulomb gauge, this agreement is already not exact at next-to-leading order Burnier:2009bk (). Calculating the potential to higher order in perturbation theory is a difficult task Brambilla:2010vq () and given the size of the strong coupling and the infrared problems in gauge theories, it is evident that non-perturbative methods within QCD, such as lattice simulations are required. The main difficulty we face is that numerical calculations are performed in imaginary time without direct access to dynamical quantities, such as .

In Ref. Rothkopf:2011db () a strategy was laid out how to evaluate the real-time definition Eq. (1) using Euclidean lattice QCD simulations. It is based on a spectral decomposition

where denotes the Euclidean time Wilson loop accessible on the lattice. The above can be combined with Eq.(1) to yield


in turn relating the values of the potential to the spectral function , which can in principle be obtained from lattice QCD.

The first practical challenge lies in obtaining the function in Eq. (2) from a finite lattice QCD dataset with statistical errors. Extracting from it continuous spectral features is an inherently ill-posed problem, which however can be given meaning by the use of Bayesian inference. In this well established statistical approach, additional prior information is used to select a unique solution from an otherwise undetermined fit. Unfortunately the standard methods, such as the Maximum Entropy Method (MEM) or extended MEM Asakawa:2000tr () have been shown Burnier:2013fca () to be unreliable in reproducing the very narrow peak structures present in the spectrum. Only recently it has become possible to faithfully reconstruct the functional form of such spectral functions from Euclidean correlator data with the advent of an improved Bayesian approach detailed in Ref. Burnier:2013nla (). In particular it is devoid of the convergence problems inherent in the MEM.

The second challenge is related to the Fourier transform and the infinite time limit in Eq. (2). Ref. Rothkopf:2011db () noted that the lowest lying spectral peak will dominate the late-time evolution and hence represents the potential contributions. A thorough understanding of the timescales influencing the shape of this spectral feature was however obtained only later in Ref. Burnier:2012az (). It was established on general grounds that if a potential picture is applicable, we have to expect a skewed Lorentzian as the lowest lying peak in the spectrum

Its position and width correspond to the values of the real- and imaginary part respectively, which can be obtained from a fit to the reconstructed spectra (i.e. with and as parameters). The practicability of this strategy has been verified within HTL perturbation theory in Burnier:2013fca () and it was found that for the purpose of extracting the potential, the Wilson loop can possibly be replaced by other observables such as the Wilson line correlators in Coulomb gauge. One particular benefit of this observable is the absence of cusp divergences and a much improved signal to noise ratio on the lattice.

Figure 1: Spectral reconstruction: On-axis Wilson line correlator data (top) at and (bottom) the spectral functions obtained by the new Bayesian reconstruction method.

Equipped with these technical and conceptual improvements we proceed to extract the temperature dependence of the static in-medium interquark potential in a purely gluonic as well as for the first time in a full QCD medium. We generated quenched QCD configurations, based on the naive anisotropic Wilson action in a fixed scale approach, i.e temperature is changed between by modifying the number of temporal lattice points (see Tab. 1). Our choice of corresponds to a relatively fine lattice spacing of fm Asakawa:2003re (), which together with a spatial extend of allows us to access both the Coulombic part of the potential, as well as those distances at which it is already screened at . We use an anisotropy of , since for a reliable determination of , a large number of points in temporal direction is required.

After fixing our configurations to Coulomb gauge, we measure the Wilson line correlators (see Tab. 1 for the number of measurements ). Their values are obtained along each spatial axis (see e.g. top of Fig. 1), on the square- and cubic diagonals. Part of the finite lattice spacing artifacts are removed by correcting the spatial distances from a comparison of the free lattice propagator with the continuum Luscher:1995zz ().

SU(3): 24 32 40 48 56 64 72 80 96
[MeV] 839 629 503 419 360 315 280 252 210
3270 2030 1940 1110 1410 1520 860 1190 1800
QGP: 6.8 6.9 7 7.125 7.25 7.3 7.48
[MeV] 148 164 182 205 232 243 286
a [fm] 0.111 0.1 0.09 0.08 0.071 0.068 0.057
1295 1340 1015 840 1220 1150 1130
Table 1: Lattice QCD configurations: (top) quenched SU(3) on anisotropic lattices with fm and MeV. (bottom) Isotropic HotQCD lattices with asqtad action (MeV).

We perform the Bayesian reconstruction of the Wilson line spectra at different temperatures excluding the first and last correlator data point at to avoid overlap divergences Berwein:2012mw (). In order to not introduce a bias for the functional form, we work with a flat default model . Frequencies are discretized at along GeV with a high resolution subinterval around the lowest lying peak. This choice is large enough for the spectra to settle parallel to the default model at large and we have checked that further extending the range or the number of points does not change the outcome. A unique global solution is found based on an LBFGS minimizer with bit precision arithmetic and a step size stopping criterion of . Several of the reconstructed spectra for are shown in Fig. 1.

In the top panel of Fig. 2 the results for the real part from the position of the lowest lying spectral peak are given by colored open symbols. They are contrasted to the color singlet free energies in Coulomb gauge , obtained on the same lattices (filled gray circles). Since the raw values fall on top of each other at small distances we have shifted them for better readability. The error bars shown are obtained from the Jackknife variance resulting from repeating the reconstruction ten times excluding a different set of 10% of the underlying measurements each. The error bands (given for ) on the other hand denote the maximum variance obtained from changing three different quantities. One corresponds to a reduction of the number of datapoints along by four and eight, the second to changing the default model normalization (, ) or functional form () and the third to the reduction in signal to noise ratio by excluding 10%,20% or 30% of the available measurements. Note that because the spectral reconstruction takes into account all datapoints along , our results for are much more robust than the free energies, that rely on a single data point. On the other hand the Bayesian reconstruction suffers from a diminishing number of datapoints at increasing temperature, as seen in the errorbands.

Our main observation is that even though the data point is excluded from the reconstruction, the values of obtained at all temperatures lie close to the color singlet free energies. While the lowest temperature shows no or very weak deviation from a linearly rising potential, the values above show clear signs of Debye screening with increasing temperature. At fm we find little temperature dependence, as expected.

Figure 2: Gluonic medium: (top) The shifted real part of the static inter-quark potential (open symbols) compared to the color singlet free energies (gray circles). Error bars represent statistical uncertainty, error bands include also systematics (see main text). (bottom) (symbols) shifted and compared to the HTL predictions (solid lines).

The extraction of the imaginary part from Bayesian spectra poses an even more formidable challenge than . Its presence can be qualitatively inferred already from the Euclidean correlator (see Fig. 1 top panel), where at intermediate values a deviation from the exponential decay and a finite curvature emerges. To obtain quantitative results, the reconstruction of the lowest lying peak needs to capture not only the width, which encodes but also the overall skewed Lorentzian shape related to non-potential effects.

The novel Bayesian approach for the first time allows us to extract this functional form (see Fig. 1 bottom panel), where the MEM yielded Gaussian like features. Previous tests based on mock data from momentum regularized HTL show that to obtain values accurate to %, datasets with datapoints are required at a high precision of . If less points are available the reconstruction tends to underestimate the width, while statistical noise leads to broadening. The former effect dominates at high temperatures and at small separation distances fm where the lattice data carries small relative errors, while at larger distances the exponential suppression of the Euclidean correlator leads to an artificially broad width.

Taking these systematic effects into account, we can estimate the values of to a lie in a band which is compatible with the expectations from HTL perturbation theory at high temperature and appears to lie slightly below HTL at low temperature. The results of our first-principles investigation are consistent with the findings by a recent modeling approach based on HTL spectral functions Bazavov:2014kva ().

Next we consider the realistic setting of a thermal QCD medium containing both gluons, as well as the light and quarks. The corresponding full QCD lattices were generated by the HotQCD collaboration Bazavov:2011nk () for the study of the QCD phase structure (see Tab.1). The Bayesian reconstructions with a common are performed using steps in a numerical interval of fixed length and a high resolution interval of points to capture the lowest lying peak. Due to the high cost in generating the configurations it is currently not possible to obtain similarly large temporal extends as in quenched QCD, even with the use of supercomputers. Therefore we focus in Fig. 3 solely on the values of the real-part (colored symbols) of the potential, which are compared to the color singlet free energies (gray circles) from the same lattices. Error bars are again obtained from Jackknife variance. The error bands () result from the maximum variation among changing the number of datapoints along by one and two, changing the normalization and functional form of the default model as well as from removing 10%,20% or 30% of the statistics.

Figure 3: Quark-gluon-plasma: The real part of the static interquark potential (open symbols) compared to the color singlet free energies in Coulomb gauge (gray circles).

At temperatures below and slightly above the pseudo-critical temperature MeV on our lattices, the Bayesian reconstruction allows us to reliably determine up to physical distances of fm. The signal of the free energies at similar T is quickly lost in the much larger statistical noise. We do not observe string breaking, up to fm, most probably due to the lattice pion mass MeV still lying above the physical value. The presence of fermionic d.o.f. significantly changes the character and location of the phase transition. Debye screening is already pronounced at the MeV. Fig. 3 shows that just as in the quenched case the real part lies close to the singlet free energy.

We have measured the static in-medium inter-quark potential defined from first principles in lattice QCD simulations at phenomenologically relevant temperatures around the deconfinement transition. The real part of the complex potential is found to lie close to the color singlet free energies in Coulomb gauge and displays the expected Debye screened behavior above the phase transition. Consequently it disagrees with the internal energies, another observable that was used as potential in phenomenological studies. It would be hence interesting to use this first principles potential as input to heavy quarkonium spectral function and phenomenological studies Burnier:2007qm (). Our estimate of the imaginary part of the potential in the quenched case also shows reasonable agreement with HTL above . At it is systematically smaller, as can be expected from the Boltzmann suppressed density of pions and glue balls present there.

While the spacing of our quenched lattices is already fine, a true continuum extrapolation needs to be pursued in the future. Some conceptual hurdles exist. It is long known that the Wilson loop does not possess a continuum limit due to cusp divergences Korchemsky:1987wg (). The Wilson line correlators contain end-point divergences, which accidentally vanish in Coulomb gauge up to two loops Aoyama:1981ev () but their non-perturbative behavior towards is unknown.

A more accurate reconstruction of will also require lattices of significantly higher temporal resolution, which is currently challenging for quenched and impractical for full QCD. A concerted effort towards the tuning of dynamical QCD lattices with anisotropy would thus be greatly beneficial to this field of research.

The authors thank H. B. Meyer, M. P. Lombardo, P. Petreczky and J.-I. Skullerud for fruitful discussions. YB is supported by SNF grant PZ00P2-142524 and AR acknowledges partial support by the DFG.


  • (1) N. Brambilla et al., Eur. Phys. J. C 71, 1534 (2011).
  • (2) A. Gray, et al., Phys. Rev. D 72, 094507 (2005); G. S. Bali, S. Collins and C. Ehmann, Phys. Rev. D 84, 094506 (2011); L. Liu et al. JHEP 1207, 126 (2012); N. Brambilla et al., arXiv:1404.3723.
  • (3) E. V. Shuryak, Phys. Lett. B 78, 150 (1978); T. Matsui and H. Satz, Phys. Lett. B 178 (1986) 416; T. Hashimoto, K. Hirose, T. Kanki and O. Miyamura, Phys. Rev. Lett. 57, 2123 (1986).
  • (4) M. Asakawa and T. Hatsuda, Phys. Rev. Lett. 92, 012001 (2004).
  • (5) A. Jakovac, P. Petreczky, K. Petrov and A. Velytsky, Phys. Rev. D 75, 014506 (2007); G. Aarts, C. Allton, M. B. Oktay, M. Peardon and J. -I. Skullerud, Phys. Rev. D 76, 094513 (2007); H. Ohno et al., Phys. Rev. D 84, 094504 (2011); G. Aarts et al., JHEP 1111, 103 (2011); H. T. Ding, A. Francis, O. Kaczmarek, F. Karsch, H. Satz and W. Soeldner, Phys. Rev. D 86, 014509 (2012). G. Aarts, et al., JHEP 1407, 097 (2014); S. Kim, P. Petreczky and A. Rothkopf, arXiv:1409.3630.
  • (6) A. Adare et al. [PHENIX Collaboration], Phys. Rev. Lett. 98 (2007) 232301; Z. Tang [STAR Collaboration], J. Phys. G 38 (2011) 124107; S. Chatrchyan et al. [CMS Collaboration], JHEP 1205 (2012) 063, Phys. Rev. Lett. 107 (2011) 052302; B. Abelev et al. [ALICE Collaboration], Phys. Rev. Lett. 109 (2012) 072301.
  • (7) S. Nadkarni, Phys. Rev. D 34, 3904 (1986); C. Y. Wong, Phys. Rev. C 72, 034906 (2005); O. Kaczmarek, PoS CPOD 07, 043 (2007); H. Satz, J. Phys. G 36, 064011 (2009).
  • (8) J. L. Albacete, Y. V. Kovchegov and A. Taliotis, Phys. Rev. D 78, (2008) 115007; T. Hayata, K. Nawa and T. Hatsuda, Phys. Rev. D 87, 101901 (2013).
  • (9) Y. Koma, M. Koma and H. Wittig, Phys. Rev. Lett. 97, 122003 (2006); T. Kawanai and S. Sasaki, Phys. Rev. Lett. 107, 091601 (2011).
  • (10) N. Brambilla, A. Pineda, J. Soto and A. Vairo, Rev. Mod. Phys. 77, 1423 (2005).
  • (11) M. Laine, O. Philipsen and M. Tassler, JHEP 0709, 066 (2007); A. Beraudo, J. -P. Blaizot, C. Ratti, Nucl. Phys. A 806, 312 (2008); N. Brambilla, J. Ghiglieri, A. Vairo, P. Petreczky, Phys. Rev. D 78, 014017 (2008).
  • (12) Y. Burnier, M. Laine and M. Vepsalainen, JHEP 1001, 054 (2010).
  • (13) N. Brambilla, M. A. Escobedo, J. Ghiglieri and A. Vairo, JHEP 1112 (2011) 116; JHEP 1305 (2013) 130.
  • (14) A. Rothkopf, T. Hatsuda and S. Sasaki, Phys. Rev. Lett. 108 (2012) 162001.
  • (15) M. Asakawa, T. Hatsuda and Y. Nakahara, Prog. Part. Nucl. Phys. 46 (2001) 459; A. Rothkopf, J. Comput. Phys. 238 (2013) 106.
  • (16) Y. Burnier and A. Rothkopf, Phys. Rev. D 87 (2013) 11, 114019; PoS LATTICE 2013 (2013) 491.
  • (17) Y. Burnier and A. Rothkopf, Phys. Rev. Lett. 111 (2013) 18, 182003; PoS LATTICE 2013 (2013) 490.
  • (18) Y. Burnier and A. Rothkopf, Phys. Rev. D 86 (2012) 051503.
  • (19) S. Necco and R. Sommer, Nucl. Phys. B 622 (2002) 328.
  • (20) M. Berwein, N. Brambilla, J. Ghiglieri and A. Vairo, JHEP 1303 (2013) 069.
  • (21) A. Bazavov, Y. Burnier and P. Petreczky, arXiv:1404.4267 [hep-lat].
  • (22) A. Bazavov et al., Phys. Rev. D 85, 054503 (2012)
  • (23) Y. Burnier, M. Laine, M. Vepsalainen, JHEP 0801 (2008) 043; A. Mocsy and P. Petreczky, Phys. Rev. D 77 (2008) 014501 A. Mocsy, P. Petreczky and M. Strickland, Int. J. Mod. Phys. A 28 (2013) 1340012; M. Strickland, Int. J. of Mod. Phys. A, Vol. 28, 1340012 (2013). A. Rothkopf, JHEP 1404, 085 (2014)
  • (24) G. P. Korchemsky and A. V. Radyushkin, Nucl. Phys. B 283 (1987) 342; R. A. Brandt, A. Gocksch, M. A. Sato and F. Neri, Phys. Rev. D 26 (1982) 3611;
  • (25) S. Aoyama, Nucl. Phys. B 194 (1982) 513.
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description