Heavy quark momentum diffusion from the lattice using gradient flow
Luis Altenkort, Alexander M. Eller, Olaf Kaczmarek, Lukas Mazur, Guy D. Moore, Hai-Tao Shu
HHeavy quark momentum diffusion from the lattice using gradient flow
Luis Altenkort , Alexander M. Eller , O. Kaczmarek , , Lukas Mazur , Guy D. Moore , H.-T. Shu Fakult¨at f¨ur Physik, Universit¨at Bielefeld, D-33615 Bielefeld, Germany Institut f¨ur Kernphysik, Technische Universit¨at DarmstadtSchlossgartenstraße 2, D-64289 Darmstadt, Germany Key Laboratory of Quark and Lepton Physics (MOE) and Institute of Particle Physics,Central China Normal University, Wuhan 430079, China (Dated: September 30, 2020)We apply the gradient flow on a color-electric two-point function that encodes the heavy quarkmomentum diffusion coefficient. The simulations are done on fine isotropic lattices in the quenchedapproximation at 1 . T c . The continuum extrapolation is performed at fixed flow time followed by asecond extrapolation to zero flow time. Perturbative calculations of this correlation function underWilson flow are used to enhance the extrapolations of the non-perturbative lattice correlator. Thefinal estimate for the continuum correlator at zero flow time largely agrees with one obtained froma previous study using the multi-level algorithm. We perform a spectral reconstruction based onperturbative model fits to estimate the heavy quark momentum diffusion coefficient. The approachwe present here yields high-precision data for the correlator and is also applicable for actions withdynamical fermions. I. INTRODUCTION
The Yang-Mills gradient flow has proven to be a power-ful tool for gauge theories since it was proposed in [1–5].Its applications on the lattice can be useful in severalways, for example reducing high frequency backgroundnoise from gauge configurations or setting the physicalscale in lattice QCD [6–8]. Observables calculated undergradient flow are claimed to be automatically renormal-ized at sufficiently large flow time if the continuum limitis taken [3, 9]. This is a convenient property when deal-ing with operators whose renormalizations are trouble-some, for instance the energy momentum tensor, whichhas already been studied under flow [10–13]. Even moreapplications of the gradient flow can be found in [14].Studies on transport coefficients like the heavy quarkmomentum diffusion coefficient, whose knowledge is ofphenomenological interest, also benefit from the gradi-ent flow. Non-perturbative determinations of such coef-ficients require high precision estimates of correspondingtwo-point correlation functions whose signals are oftenovershadowed by high frequency background noise evenat large Monte Carlo sample sizes. Techniques like themulti-level algorithm [15] and link-integration [16] canameliorate the problem but are not applicable for ac-tions with dynamical quarks. For now the gradient flowseems to be the only way to obtain high precision dataof correlation functions in simulations with dynamicalquarks. Discussions on how to flow dynamical quarksand their applications in dynamical QCD can be foundin [5, 13, 17].In this paper we want to study how the well-knowncorrelation function of color-electric fields which containsinformation about heavy quark transport behaves underYang-Mills gradient flow and how to perform its con-tinuum and flow-time-to-zero extrapolation. This color-electric correlator [18] has been intensively studied in per-turbation theory [18–21] and also non-perturbatively [22– 25]. For this study we restrict ourselves to pure SU(3)gauge theory at a temperature of T ≈ . T c . We com-pare our results with a previous non-perturbative mea-surement of the correlator from [24] that was obtained inthe same setting but with the multi-level algorithm as anoise reduction method. This allows us to cross-check theresults obtained from the gradient flow approach. Theextrapolated correlator in the continuum limit at zeroflow time is then used for spectral function reconstruc-tion, for which we use multiple theoretically motivatedfitting models, similar to what was done in [24]. Finallywe compare the heavy quark momentum diffusion coeffi-cient obtained from the extracted spectral function withthose from other studies [24, 26].The paper is structured as follows: we start by recallingthe definition of the gradient flow and how it is realizedon the lattice (section II). Afterwards we briefly presentsome perturbative calculations of the color-electric corre-lator at non-zero flow time and how we can use them toenhance the non-perturbative lattice data (Section III).Section IV is then devoted to the analysis of this dataand how we carry out the double extrapolation. The fi-nal estimate for the renormalized continuum correlatorat zero flow time is then used to obtain the heavy quarkmomentum diffusion coefficient through spectral recon-struction in section V. We summarize our findings insection VI. II. GRADIENT FLOW
The Yang-Mills gradient flow evolves the gauge field tosome flow time τ F along the gradient of the gauge action.The flowed gauge field B µ ( τ F , x ) is defined through B µ ∣ τ F = = A µ , ∂B µ ∂τ F = D ν G νµ , (1) a r X i v : . [ h e p - l a t ] S e p where A µ is the ordinary gauge field. In order to ful-fill gauge invariance one can easily construct the flowedcovariant derivative and field strength following the stan-dard Yang-Mills gauge theory: D µ = ∂ µ + [ B µ , ⋅ ] ,G µν = ∂ µ B ν − ∂ ν B µ + [ B µ , B ν ] . (2)In perturbation theory, a composite local operator O ( x, τ F ) consisting of gauge fields can be expanded in τ F as a superposition of the renormalized operator definedat zero flow time [4]: O ( x, τ F ) ———→ τ F → ∑ i c i ( τ F ) O Ri ( x ) . (3)Here c i are some coefficients that can be calculated per-turbatively. In order to obtain the renormalized operatorin the original theory one needs to invert Eq. 3 and goback to zero flow time. We will show how to carry outthis procedure in practice in section IV D.To leading order a perturbative calculation in D di-mensions [3] provides a simple explanation for the effectof the flow on the gauge field: B µ ( x, τ F ) = ∫ d D y K τ F ( x − y ) A µ ( y ) ,K τ F ( z ) = ∫ d D p ( π ) D e ipz e − τ F p = e − z / τ F ( πτ F ) D / . (4)The flow averages the gauge field with a local Gaussiankernel. The mean-square radius of the Gaussian distri-bution is √ τ F and is sometimes called the “flow radius.”On the lattice we need to use the discretized coun-terparts of the continuum flow equations. Specifically,we use a Symanzik improved version of the flow called Zeuthen flow [27], which eliminates O( a ) cutoff effectsthat are present in the ordinary Wilson flow. We solveEq. 1 numerically using an adaptive step-size Runge-Kutta method for Lie groups [28].How much flow can be applied to a gauge field beforemeasuring a correlation function on the lattice? Considerthe correlator of two operators whose temporal separa-tion is τ . On the one hand, the flow radius should be largeenough to suppress lattice discretization effects, and onthe other hand, we need to make sure that the separationbetween operators is larger than the flow radius, so that the correlation between them is not contaminated. Thisleads to a general allowed flow time range of a ≲ √ τ F ≲ τ . (5)We can make the latter limit more rigorous by performinga perturbative determination of the flow radius where aspecific correlator starts to suffer contamination due tothe overlap of the operator smearing radii. We will dothis for the color-electric correlator below. III. THE COLOR-ELECTRIC CORRELATOR INPERTURBATION THEORY
The correlator that we want to study under flow wasfirst proposed in [18] and is related to the thermalizationrate of a heavy quark in a hot plasma. It arises from theinfinite quark mass limit of the meson current-currentcorrelator. The advantage of studying this specific cor-relator is that the spectral function encoded in it hasno sharp transport peak, in contrast to the perturbativebehavior for current-current and stress-stress correlators.As a result, the relevant transport coefficient is much eas-ier to obtain. The correlator is a product of electric fieldsthat sit on a Polyakov loop, which is why it is referred toas the color-electric correlator. It is defined as [18] G ( τ ) = − ∑ i = ⟨ Re Tr [ U ( β, τ ) E i ( τ, ⃗ ) U ( τ, ) E i ( , ⃗ )]⟩ ⟨ Re Tr [ U ( β, )]⟩ , (6)where β = / T is the temporal extent which is equal tothe inverse temperature, and U ( τ , τ ) is a Wilson line inthe Euclidean time direction. E i is the color-electric fieldin the geometrical normalization, which we discretize fol-lowing [18]: E i ( τ, ⃗ x ) = U i ( τ, ⃗ x ) U ( τ, ⃗ x + ˆ i )− U ( τ, ⃗ x ) U i ( τ, ⃗ x + ˆ4 ) . (7)Here, U i ( τ, ⃗ x ) is an SU(3) link variable in the i th direc-tion.The leading-order perturbative calculation in contin-uum of the correlator under gradient flow was performedin [20]. We perform the analog calculation on the latticehere in order to illustrate the behavior of the correlatorunder Wilson flow compared to its continuum counter-part.The LO perturbative lattice correlator under Wilson flow reads G normlatt ( τ, τ F ) ≡ G LOlatt ( τ, τ F ) g C F = − a β β − ∑ n =− β cos ( πnτ T ) (8) × ⎡⎢⎢⎢⎢⎢⎣ e − τ F T N τ sin ( π nNτ ) e − τ F T N τ ( I ( τ F T N τ )) − ∞ ∫ τ F T N τ dx e − x sin ( π nNτ ) e − x ( I ( x )) { I ( x ) − I ( x )} ⎤⎥⎥⎥⎥⎥⎦ , . . . . . . τ T G norm τ F ( τ ) √ τ F T =0.00, cont0.00, latt0.05, cont0.05, latt0.10, cont0.10, latt FIG. 1. The color-electric correlator at different flow timescalculated perturbatively to leading order in the continuum[20] and on the lattice (Eq. 8) with N τ =
20. The vertical linesindicate the point where the flowed and unflowed (continuum)results begin to agree to better than 1%. For the lattice cor-relators this distance needs to be shifted to the right by aboutone lattice spacing. where I , I are modified Bessel functions of the firstkind and C F = ( N c − )/ N c = / g C F and equipit with the superscript “norm” since we will later use itto normalize our non-perturbative lattice measurements.In Figure 1 we show an example for different flow timesin the continuum and on a lattice with N τ =
20. Thevertical dashed line in this figure is the lower limit ofthe distance for which the flowed continuum correlatordeviates by less than 1% from its non-flowed counterpartand is given by τ ≳ √ τ F [20]. For the lattice versionof this limit one needs to shift the allowed distance byabout one lattice spacing, and together with Eq. 5 oneobtains a serviceable flow interval of a ≲ √ τ F ≲ τ − a . (9)Because of this limitation one can only hope to extractthe large-distance behavior via the gradient flow method.Fortunately this is also the region that carries the infor-mation about the heavy quark momentum diffusion co-efficient.From Figure 1 we can see that the correlator growsas τ − at small separations, and therefore spans severalorders of magnitude in value. In order to improve vis-ibility we can remove this dominant behavior from thenon-perturbative lattice data that we present later bynormalizing it with the perturbative leading-order cor-relator. The normalized data is also convenient for theextrapolations as its values are all of the same order ofmagnitude. Correlators calculated on the lattice suffer from cutoffeffects. One can remove the leading-order contributionof these effects from the non-perturbative lattice databy utilizing the perturbative lattice and continuum cor-relators (see Eq. 8 and [20]). This technique is knownas tree-level improvement [29]. In appendix A 2 we ex-plain in detail how to perform the improvement for thecolor-electric correlator under flow. Note that in ourstudy the non-perturbative correlator is obtained usingZeuthen flow, while the perturbative lattice calculationsare done with Wilson flow. The main point of appendixA 2 is that it is reasonable to improve the Zeuthen-flowedcorrelator by using the non-flowed leading-order latticecorrelator. Since we want to normalize the data with theleading-order continuum correlator, we end up with thedimensionless ratio (cf. Eq. A6 and A7) G latt τ F ( τ ) G normlatt τ F = ( τ ) , (10)on which the continuum and flow-time-to-zero extrapo-lations will be performed. IV. FLOWED COLOR-ELECTRICCORRELATOR ON THE LATTICE
In this section we present our measurements and analy-sis of the color-electric correlator on fine isotropic latticesat T ≈ . T c . We perform the continuum extrapolationat fixed flow time and subsequently the extrapolation to τ F = A. Lattice setup
A summary of the parameters of the gauge configura-tions used in this study can be found in Table I. All con-figurations are generated in the quenched approximationusing heatbath and overrelaxation updates. The gaugeaction is the standard Wilson action. After thermaliza-tion (5000 heatbath sweeps), we save a configuration af-ter every 500 combined sweeps, where one such sweepconsists of one heatbath and four overrelaxation sweeps.The lattice boundary conditions are periodic for all di-rections.In order to compare our correlator measurements withprevious results from [24], we use the same lattice di-mensions and β -values to set the temperature to about1 . T c . The lattice spacing a is determined via the Som-mer scale r [29] with parameters taken from [30] andupdated coefficients from [31]. We use the state-of-the-art value r T c = . ( ) [30].For statistical mean and error estimation of observ-ables measured on the configurations we perform a boot-strap analysis with 10000 samples for each lattice. cont. N τ = 36 N τ = 30 N τ = 24 N τ = 200 . . . . τ T . . . . . . . . . G latt τ F ( τ ) G normlatt τ F =0 ( τ ) √ τ F T = 0.060 cont. N τ = 36 N τ = 30 N τ = 24 N τ = 200 . . . . τ T . . . . . . . . . G latt τ F ( τ ) G normlatt τ F =0 ( τ ) √ τ F T = 0.100 FIG. 2.
Left : The color-electric correlator at T ≈ . T c on different lattices at a small flow time ( √ τ F T = . a = /( N τ T ) with N τ = Right : The same but at √ τ F T = . B. Cutoff effects and signal-to-noise ratio underflow
Calculations of the color-electric correlator withoutgradient flow exhibit huge statistical errors even for thelarge amounts of gauge configurations used in this study(see Table I). Beyond the first few lattice spacings theobtained data is almost pure noise at zero flow time. Byincreasing the flow time, high-frequency fluctuations ofthe gauge fields are suppressed and the signal-to-noiseratio of the correlator is expected to increase. Addition-ally, cutoff effects are expected to decrease with increas-ing flow time.The minimum amount of flow that is needed to pro-duce renormalized operators is √ τ F T ≈ aT = N − τ . ViaEq. 9 one can convert this minimum flow time into a lowerbound for the separation τ T . Because of this we hence-forth exclude the N τ =
16 lattice (see Table I) from theanalysis. The limits are then dictated by the N τ =
20 lat- a (fm) a − (GeV) N σ N τ β T / T c β values, tem-perature and number of configurations generated for thiswork. The lattice spacing a is determined via the Sommerscale r [29] with parameters taken from [30] and updatedcoefficients from [31]. We use r T c = . ( ) [30]. Thecoarsest lattice ( N τ =
16) does not enter in the continuumextrapolation. tice, and so we have to flow up to at least √ τ F T = . τ T ≥ .
2. These limits are also applied to the finer latticesand will carry over to the continuum and flow-time-to-zero extrapolations.In Figure 2 we show the tree-level improved correlatorson all lattices at one intermediate ( √ τ F T = .
06) andone larger flow time ( √ τ F T = . N τ = τ F T . . . . . . . G latt τ ( τ F ) G normlatt τ,τ F =0 τ T =0.5000.4720.4440.4170.3890.3610.3330.3060.278 0.2500.2220.1940.1670.1390.1110.0830.0560.0280.0 0.05 0.08 0.10 0.12 0.13 0.14 0.15 √ τ F T FIG. 3. The color-electric correlator as a function of flowtime on the finest available lattice ( N τ =
36) at T ≈ . T c .The grey error bands indicate the underlying data; to reduceclutter they are hidden if statistical errors would be largerthan 5%. In retrospect we only consider the flow times forthe flow extrapolation where the data points are shown ex-plicitly (see section IV D). The dashed vertical line is theminimum amount of flow that is needed for the coarsest con-sidered lattice ( N τ =
20, see Table I). The horizontal positionsof the various markers depict the flow limits from Eq. 9 (with a = /( N τ T ) , N τ = plete and the correlation then behaves almost linearly forsome time, before it is then substantially contaminated.After a continuum extrapolation we will exploit the seem-ingly linear region to extrapolate the correlator to zeroflow time. The dominant flow-time-dependent contribu-tions in this region will essentially consist of graduallyaccruing contact terms of the operator. In principle theseexist as soon as the flow time is non-zero, and their con-tribution is enhanced by the fact that the operator is nottruly local due to the discretization of the electric fields(Eq. 7). C. Continuum extrapolation
The physical correlation function is at vanishing latticespacing and flow time, so we need to perform a doubleextrapolation. But it is important to perform the extrap-olations in the right order: the continuum extrapolationshould be performed first, and the flow-time-to-zero ex-trapolation should be carried out on the continuum corre-lation functions. Otherwise, lattice-spacing issues whichare ameliorated by flow re-emerge in the flow-time-to-zero extrapolation.The continuum extrapolations should only be per-formed at the separations of the finest lattice that isavailable, as the information is most reliable there. Alldata from coarser lattices need to be interpolated to these . . . . . . N − τ . . . . . . G latt τ,τ F G normlatt τ,τ F =0 τ T =0.5000.4720.444 0.4170.3890.361 0.3330.306 √ τ F T = 0 . FIG. 4. Continuum extrapolation of the color-electric corre-lator at √ τ F T = .
085 using an
Ansatz linear in N − τ . Thelattices have temporal extents N τ ∈ { , , , } . separations, for which we use cubic splines. For a givenlattice and flow time we only interpolate inside the in-terval given by Eq. 9, with the exception of the one datapoint that lies just outside of the lower separation limit.By also including this slightly contaminated point we canalways interpolate at the lower limit without explicitlyusing the contaminated data, which extends the usableflow range. At the left boundary of the spline the sec-ond derivative is set to zero (natural condition) and atthe middle point ( τ T = .
5, which in practice is the rightboundary) the first derivative is set to zero (symmetry).We perform the interpolations on every bootstrap sampleand take the bootstrap mean and standard deviation ofthe interpolated points for statistical error estimation.Once the interpolated values are obtained, the contin-uum extrapolation is carried out, independently for eachdistance τ and flow time τ F , by performing a weightedfit on the data with the Ansatz G latt τ,τ F ( N τ ) G normlatt τ,τ F = ( N τ ) = m ⋅ N − τ + b (11)with fit parameters m and b = G cont τ,τ F / G normcont τ,τ F = . The Wilsonaction that we use to generate the gauge configurationshas leading discretization errors of order a , which ex-plains this choice. We estimate the statistical error byperforming the weighted fit on substitute bootstrap sam-ples that we draw from a Gaussian around the bootstrapmean of the interpolations with their error as the width.The bootstrap mean and standard deviation of the sub-stitute samples then serve as the statistical estimates. InFigure 4 we show the continuum extrapolation for oneintermediate flow time ( √ τ F T = . τ T =0.5000.4720.4440.4170.3890.361 0.3330.3060.2780.2500.2220 . . . . . τ F T . . . . . . . G cont τ ( τ F ) G normcont τ,τ F =0 √ τ F T FIG. 5. Flow-time-to-zero extrapolation of the continuumcolor-electric correlator at fixed τ T using a linear
Ansatz . Therange of flow times at each separation is restricted by Eq. 9and the statistical precision (see section IV D).
D. Flow-time-to-zero extrapolation
The continuum correlators at fixed flow time still needto be extrapolated back to τ F =
0. Because of the limitedknowledge we have of the functional form of the floweffect we use the simple linear
Ansatz G cont τ ( τ F ) G normcont τ,τ F = = m ⋅ τ F + b, (12)where m and b = G cont τ / G normcont τ are tau -dependent fit pa-rameters. For statistical error estimation we proceed asin the continuum extrapolation.Based on the continuum correlators we estimate byhand for which separations and flow times this Ansatz isreasonable (see Figure 5 and also Figure 3):For the lower boundary of the flow time the criterion isthe statistical precision of the data. There are two waysto improve the precision: 1. apply more flow, 2. increasethe number of gauge configurations. For the larger dis-tances where one needs more configurations for the samestatistical accuracy (in comparison to the smaller dis-tances), the minimum amount of flow from Eq. 5 is notenough to give sufficiently small statistical errors on thedata. Therefore we only use continuum data points in theflow-time-to-zero extrapolation that are flowed enoughso that they have a maximum relative error of 1% ⋅ τ T .This value is chosen to yield a lower boundary for theflow time where a linear Ansatz is reasonable. This alsoimplies that the minimal separation for which we canobtain renormalized correlator data is shifted. For theupper boundary of the flow time we are limited by Eq. 9and by the number of data points of the coarsest lattice; ifthere are less than two non-contaminated data points forany given lattice, we cannot perform the interpolation
Gradient flow methodMulti-level method0.10 0.20 0.30 0.40 0.5 τ T . . . . . G cont ( τ ) G normcont ( τ ) FIG. 6. Non-perturbatively renormalized continuum color-electric correlator at zero flow time obtained from the gradientflow method (with lattice setup from Table I) in comparisonwith revised continuum correlator from multi-level method(perturbatively renormalized to NLO [32]) using lattice datafrom [24]. and continuum extrapolation. For the coarsest latticethis turns out to be at a flow time of √ τ F T > . τ T ∈ [ / , . ] . In Figure 6 we show the final renormal-ized continuum correlator at zero flow time. The errorsshown here only depict statistical uncertainties. We alsocarry out a new continuum extrapolation of the latticecorrelators from a previous study [24] which used themulti-level algorithm and link integration techniques toobtain a signal for the color-electric correlator. The oldmethod only works in pure gauge theory and needs tobe renormalized approximately with the help of pertur-bation theory [32]. Note that the statistics used in [24]are much smaller than the ones in this study (Table I),especially for the finer lattices.The two results agree in their overall shape, while theslight offset can be explained by missing higher-order con-tributions to the renormalization of the multi-level result,the systematic uncertainty introduced by the flow extrap-olation, and the difference in lattices and statistics.A note on systematic uncertainties: In principle, thecorrelator data shown in Figure 6 is subject to system-atic uncertainties through scale setting, auto-correlationbetween gauge configurations and cross-correlation be-tween the correlator distances, the non-flowed tree-levelimprovement of the Zeuthen-flowed correlators, and fi-nally the interpolation and extrapolation Ans¨atze . V. HEAVY QUARK MOMENTUM DIFFUSIONCOEFFICIENTA. Spectral function of the color-electric correlator
It was shown in [18] that the heavy quark momentumdiffusion coefficient can be obtained from κ = lim ω → Tω ρ ( ω ) , (13)where ρ ( ω ) is the spectral function of the color-electriccorrelator, which is related to the Euclidean correlatorwe compute here through the integral relation G ( τ ) = ∫ ∞ d ωπ ρ ( ω ) cosh ( ω ( τ − T )) sinh ( ω T ) . (14)In the following we will calculate the leading-order con-tinuum spectral function under gradient flow, which willlead to a discussion about the Kramers-Kronig relationand flowed retarded correlators.The (continuum, zero flow time) leading-order spectralfunction is well known [21]. But we did not find a study ofthe spectral function for the Euclidean function at finiteflow depth, so we will consider this problem here. Westart from the leading-order correlation function at flowtime τ F , G LO ( τ, τ F ) = − g C F ⨋ K e ik n τ e − τ F K ( D − ) k n + k K = − g C F d ( π ) d / T ∑ k n e ik n τ ( k n ) d / × ( Γ ( − d / , τ F k n ) +
12 Γ (− d / , τ F k n )) , (15)where Γ ( s, x ) = ∫ ∞ x d t t s − e − t is the incomplete Gammafunction and d the number of space dimensions. Fouriertransforming the correlator to frequency space using the(periodic) Kronecker delta function ∫ β d τ e ik n τ = βδ k n , leads to˜ G LO ( ω n , τ F ) = − g C F d ( π ) d / ( ω n ) d / × ( Γ ( − d / , τ F ω n ) +
12 Γ (− d / , τ F ω n )) . (16)Setting the dimension to its physical value d = ρ ( ω ) = Im [ lim (cid:15) → + ˜ G ( ω n → − iω + (cid:15) )] , (17) we find the spectral function to be ρ LO ( ω, τ F ) = − g C F ω ( π ) / (18) × lim (cid:15) → + Re [ Γ ( − / , τ F (− iω + (cid:15) ) )+
12 Γ ( − / , τ F (− iω + (cid:15) ) )] . The real parts of the incomplete Gamma functions do notvary with their incomplete argument. The value of thereal part of the incomplete Gamma function is thereforethe same as the ordinary Gamma function. The leadingorder spectral function is ρ LO ( ω, τ F ) = g C F π ω , (19)which is the the same result as the non-flowed spectralfunction calculated in [21]. Therefore, at leading per-turbative order, the spectral function for the flowed andfor the non-flowed Euclidean correlation function are thesame. At first this appears remarkable; two distinct Eu-clidean functions reconstruct to the same spectral func-tion, which implies that Eq. 14 must not apply to theflowed correlation function.The connection between the retarded correlator andthe spectral function is made by the Kramers-Kronig re-lation (see [33]). The real and imaginary parts of a mero-morphic function, which is analytic in the closed upperhalf-plane, are related via Kramers-Kronig integral for-mulas if the function vanishes fast enough for a growingcomplex argument. The mathematical requirements ofthe function can be translated to the necessity of causal-ity of the underlying physical theory.The introduction of gradient flow, that is, replacing lo-cal operators with smeared-out relatives, breaks energy-momentum conservation through new contact terms.The commutator correlator of smeared-out operatorsno longer vanishes outside the light cone, thus break-ing causality. Alternatively, we can appeal to Eq. 4,where e − τ F p represents an exponential suppression forEuclidean 4-momentum p , but is an exponential enhance-ment for timelike Minkowski 4-momenta. Mathemat-ically speaking, the exponential suppression along theEuclidean time axis τ will turn to exponential growthafter Wick rotation. This happens during the analyticcontinuation of the Euclidean to the retarded correlator,which no longer fulfills the requirements needed for theKramers-Kronig relation, thereby breaking its connectionto the spectral function.This discussion proves that spectral reconstruction canonly be attempted on data which has been extrapolatedto zero flow time (as well as to the continuum). As wehave already emphasized, this extrapolation should beperformed after extrapolating to the continuum limit.With the extrapolated data we can apply spectral recon-struction techniques to estimate the heavy quark momen-tum diffusion coefficient κ , which we will do in the nextsection. − . − . − . − . . . . . .
25 0 . .
35 0 . .
45 0 . τTG cont − G model G normcont from ρ ( α,a )model from ρ ( α,b )model . . ω/TρωT φ ( a ) UV from ρ ( α,a )model from ρ ( α,b )model FIG. 7.
Left : The differences between the double-extrapolated correlators and the fit correlators for two models.
Right : Thefit spectral functions from different models (with n max = B. Extraction of the heavy quark momentumdiffusion coefficient
There are many approaches for a spectral reconstruc-tion [34–39]. Here we will proceed similarly to [24] andreconstruct the spectral function based on theoreticallymotivated model fits of the correlator data obtained fromthe gradient flow method that is shown in Figure 6. Thefit ansatz is a combination of two parts: one is valid in thesmall frequency regime ( ω ≪ T ) and the other is valid inthe large frequency regime ( ω ≫ T ). In the intermediateregime an interpolation is required. The small frequencypart can be described by infrared asymptotics [18], φ IR ( ω ) ≡ κω T , (20)and the large frequency part can be calculated in pertur-bation theory with the help of asymptotic freedom [21,40], φ ( a ) UV ( ω ) ≡ g ( ¯ µ ω ) C F ω π , ¯ µ ω ≡ max ( ω, πT ) , (21)or, when considering the higher order correction, φ ( b ) UV ( ω ) ≡ φ ( a ) UV ( ω )[ + ( r + r (cid:96) ) a s ( ¯ µ ω )] . (22) C F , r , r , (cid:96) can be found in [21, 24] and g , a s areevaluated using four-loop running [41]. Eq. 21 and Eq. 22will be treated equally in our fits as explained in [24]. Toaccount for the corrections when incorporating φ IR and φ UV , which are valid only in their own regimes, we needto introduce two trigonometric functions, e ( α ) n ( y ) ≡ sin ( πny ) , e ( β ) n ( y ) ≡ sin ( πy ) sin ( πny ) , (23)where y ≡ x + x and x ≡ ln ( + ωπT ) . Putting all thistogether yields three types of models as in [24]; however, in this work we only adopt the most theoretically justifiedone, ρ ( µ,i ) model ( ω ) ≡ [ + n max ∑ n = c n e ( µ ) n ( y )]√[ φ IR ( ω )] + [ φ ( i ) UV ( ω )] , (24)where µ ∈ { α, β } , i ∈ { a, b } . Substituting the true spectralfunction in Eq. 14 with the model spectral function, weobtain a model correlator G model with which we define aweighted sum of squares for a fit: χ ≡ ∑ τ [ G cont ( τ ) − G model ( τ ) δG cont ( τ ) ] . (25)Here δG cont ( τ ) is the statistical error of the lattice data G cont ( τ ) which we show in Figure 6.Now it is clear that the parameters to fit are κ / T and c n . As in [24], two strategies are considered. The firststrategy is just a normal fit and in the second strategywe constrain that at the maximum frequency the spectralfunction reproduces its UV asymptotics. In strategy I wedo not stop after 200 iterations as in [24] as we have seenthat this criterion makes almost no difference. For bothstrategies in our fit routine we set the absolute tolerancefor χ to 0.0001 and the maximum number of iterationsto 2000. In our fits we limit ourselves to n max ∈ { , } and we use all the data points shown in Figure 6, that is, τ T ≥ / ρ ( α,a ) model and ρ ( α,b ) model with n max = αaαbβaβb m o d e l κ/T strategy Istrategy II FIG. 8. κ / T from different fit models and strategies. Foreach model the lower data point corresponds to n max = n max = correlators are also similar. This also holds for all theother choices of models or fit strategies and we obtained χ / d.o.f. ∼ . . . . . κ / T issummarized in Figure 8. For each point the error baris the statistical uncertainty we obtain from a bootstrapanalysis and we take the spreading over various modelsas the systematic uncertainty of our estimation. Finally,based on the central values of Figure 8 we obtain a rangefor κ / T : κ / T = . . . . . , (26)which is slightly shifted to larger values compared to1 . . . . . . . . . .
64 from [26]. With theestimated range for κ / T we can also estimate the heavyquark diffusion coefficient in the heavy quark mass limit ( M ≫ πT ) according to D = T / κ [18]: DT = . . . . . . (27) VI. CONCLUSION
We calculated the color-electric correlator on fourisotropic lattices at T ≈ . T c in the quenched approxi-mation under gradient flow. We found that the gradientflow technique significantly improves the signal of thecorrelator on the lattice. It is necessary to perform botha continuum and a flow-time-to-zero extrapolation; weshow that the continuum limit must be taken first, andwe perform both extrapolations on our data within therange τ T ∈ [ / , . ] . We found that an Ansatz linearin 1 / N τ and one linear in τ F is suitable for the respec-tive extrapolations. The overall shape of the correlatorobtained via the flow method in the a → τ F → a → χ fits performed on the cor-relators using theoretically well-established models.This work affirms the validity of the gradient flow tech-nique when performing the continuum and flow-time-to-zero extrapolations and provides a method for extendingthese studies to full QCD. ACKNOWLEDGMENTS
We would like to thank Mikko Laine for useful discus-sions and previous collaboration. This work is supportedby the Deutsche Forschungsgemeinschaft (DFG, GermanResearch Foundation)-Project number 315477589-TRR211. The computations in this work were performed onthe GPU cluster at Bielefeld University. [1] R. Narayanan and H. Neuberger, JHEP , 064 (2006),arXiv:hep-th/0601210 [hep-th].[2] M. L¨uscher, Commun. Math. Phys. , 899 (2010),arXiv:0907.5491 [hep-lat].[3] M. L¨uscher, JHEP , 071 (2010), [Erratum:JHEP03,092(2014)], arXiv:1006.4518 [hep-lat].[4] M. L¨uscher and P. Weisz, JHEP , 051 (2011),arXiv:1101.0963 [hep-th].[5] M. L¨uscher, JHEP , 123 (2013), arXiv:1302.5246 [hep-lat].[6] S. Borsanyi et al. , JHEP , 010 (2012), arXiv:1203.4469[hep-lat].[7] A. Bazavov et al. (MILC), Phys. Rev. D93 , 094510(2016), arXiv:1503.02769 [hep-lat].[8] M. Dalla Brida and A. Ramos, Eur. Phys. J.
C79 , 720 (2019), arXiv:1905.05147 [hep-lat].[9] M. Luscher, PoS
LATTICE2010 , 015 (2010),arXiv:1009.5877 [hep-lat].[10] M. Kitazawa, T. Iritani, M. Asakawa, T. Hatsuda, andH. Suzuki, Phys. Rev. D , 114512 (2016).[11] M. Kitazawa, T. Iritani, M. Asakawa, and T. Hatsuda,Phys. Rev. D96 , 111502 (2017), arXiv:1708.01415 [hep-lat].[12] H. Suzuki, PTEP , 083B03 (2013), [Erratum:PTEP2015,079201(2015)], arXiv:1304.0533 [hep-lat].[13] Y. Taniguchi, S. Ejiri, R. Iwami, K. Kanaya, M. Ki-tazawa, H. Suzuki, T. Umeda, and N. Wakabayashi,Phys. Rev. D , 014509 (2017), [Erratum: Phys.Rev.D99, 059904 (2019)], arXiv:1609.01417 [hep-lat].[14] M. L¨uscher, PoS LATTICE2013 , 016 (2014), arXiv:1308.5598 [hep-lat].[15] M. L¨uscher and P. Weisz, JHEP , 010 (2001),arXiv:hep-lat/0108014 [hep-lat].[16] P. De Forcrand and C. Roiesnel, Phys. Lett. , 77(1985).[17] H. Makino and H. Suzuki, PTEP , 063B02 (2014),[Erratum: PTEP2015,079202(2015)], arXiv:1403.4772[hep-lat].[18] S. Caron-Huot, M. Laine, and G. D. Moore, JHEP ,053 (2009), arXiv:0901.1195 [hep-lat].[19] J. Casalderrey-Solana and D. Teaney, Phys. Rev. D74 ,085012 (2006), arXiv:hep-ph/0605199 [hep-ph].[20] A. M. Eller and G. D. Moore, Phys. Rev.
D97 , 114507(2018), arXiv:1802.04562 [hep-lat].[21] Y. Burnier, M. Laine, J. Langelage, and L. Mether,JHEP , 094 (2010), arXiv:1006.0867 [hep-ph].[22] H. B. Meyer, New J. Phys. , 035008 (2011),arXiv:1012.0234 [hep-lat].[23] D. Banerjee, S. Datta, R. Gavai, and P. Majumdar,Phys. Rev. D , 014510 (2012).[24] A. Francis, O. Kaczmarek, M. Laine, T. Neuhaus,and H. Ohno, Phys. Rev. D , 116003 (2015),arXiv:1508.04543 [hep-lat].[25] N. Brambilla, V. Leino, P. Petreczky, and A. vairo(TUMQCD), PoS LATTICE2019 , 108 (2019),arXiv:1912.00689 [hep-lat].[26] N. Brambilla, V. Leino, P. Petreczky, and A. Vairo,(2020), arXiv:2007.10078 [hep-lat].[27] A. Ramos and S. Sint, Eur. Phys. J. C , 15 (2016),arXiv:1508.05552 [hep-lat].[28] E. Celledoni, A. Marthinsen, and B. Owren, Future Gen-eration Computer Systems , 341 (2003), special Issueon Geometric Numerical Algorithms.[29] R. Sommer, Nucl. Phys. B , 839 (1994), arXiv:hep-lat/9310022.[30] A. Francis, O. Kaczmarek, M. Laine, T. Neuhaus,and H. Ohno, Phys. Rev. D , 096002 (2015),arXiv:1503.05652 [hep-lat].[31] Y. Burnier, H. T. Ding, O. Kaczmarek, A. L. Kruse,M. Laine, H. Ohno, and H. Sandmeyer, JHEP , 206(2017), arXiv:1709.07612 [hep-lat].[32] C. Christensen and M. Laine, Phys. Lett. B , 316(2016), arXiv:1601.01573 [hep-lat].[33] J. D. Jackson, Classical Electrodynamics (Wiley, 1998).[34] A. N. Tikhonov, A. V. Goncharsky, V. V. Stepanov,and A. G. Yagola,
Numerical methods for the solutionof ill-posed problems , Mathematics and Its Applications(Springer, Dordrecht, 1995).[35] J. Otsuki, M. Ohzeki, H. Shinaoka, and K. Yoshimi,Phys. Rev. E , 061302 (2017).[36] M. Jarrell and J. Gubernatis, Physics Reports , 133(1996).[37] Y. Burnier and A. Rothkopf, Phys. Rev. Lett. ,182003 (2013), arXiv:1307.6106 [hep-lat].[38] H.-T. Ding, O. Kaczmarek, S. Mukherjee, H. Ohno,and H. Shu, Phys. Rev. D , 094503 (2018),arXiv:1712.03341 [hep-lat].[39] G. Backus and F. Gilbert, Geophysical Journal Interna-tional , 169 (1968).[40] S. Caron-Huot, Phys. Rev. D79 , 125009 (2009),arXiv:0903.3958 [hep-ph].[41] T. van Ritbergen, J. Vermaseren, and S. Larin, Phys.Lett. B , 379 (1997), arXiv:hep-ph/9701390.[42] S. Capitani, Phys. Rept. , 113 (2003), arXiv:hep- lat/0211036 [hep-lat].[43] A. Cheng, A. Hasenfratz, Y. Liu, G. Petropoulos, andD. Schaich, JHEP , 137 (2014), arXiv:1404.0984 [hep-lat].[44] A. Hasenfratz, D. Schaich, and A. Veernala, JHEP ,143 (2015), arXiv:1410.5886 [hep-lat]. APPENDIXAppendix A: Tree-level improvement and flow1. Derivation of leading-order lattice correlatorunder Wilson flow
In this section we present the derivation of the leading-order lattice correlator under Wilson flow. The latticegluon propagator under Wilson flow has the form ⟨ B aµ ( K, τ F ) B bν ( Q, τ F )⟩ = δ ab ( π ) δ ( ) ( K + Q )× K [( δ µν − ˜ K µ ˜ K ν ˜ K ) e − τ F ˜ K + ξ ˜ K µ ˜ K ν ˜ K e − α τ F ˜ K ] , where ξ and α are the gauge and flow gauge fixing pa-rameters and we used the tilde short-hand notation˜ k µ = a sin ( ak µ ) . (A1)In the following we will use both Feynman and flowFeynman gauge.The leading-order expression reads G LOlatt ( τ, τ F ) = − g C F ⨋ K e ik n τ e − τ F ˜ K [ − ∑ i ˜ k i ˜ K ] , (A2)which is analogous to the continuum expression but ex-changing ˜ k → k . The spatial integration can be per-formed analytically π ∫ − π d k ( π ) e − τ F a ∑ i sin ( ki ) = e − τ F a ( I ( τ F a )) , π ∫ − π d k ( π ) ∑ i sin ( k i ) e − x ∑ i sin ( ki ) = e − x ( I ( x )) { I ( x ) − I ( x )} , where I and I are modified Bessel functions of thefirst kind. Inserting this into Eq. A2 and introducing aSchwinger parameter in order to rewrite the denominatorleads to the expression given in Eq. 8.
2. Flow-time-to-zero extrapolation and comparisonof tree-level improvement methods under flow
Correlators calculated on the lattice suffer from cutoffeffects. One can extract and account for the leading-order1 τ T =0.5000.438 0.3750.3120 . . . . . τ F . . . . . . . G Zeuthenlatt τ ( τ F ) G normlatt τ,τ F =0 √ τ F T τ T =0.5000.438 0.3750.3120 . . . . . τ F . . . . . . . G Wilsonlatt τ ( τ F ) G normlatt τ ( τ F ) √ τ F T FIG. 9.
Left : Naive flow-time-to-zero extrapolation of the
Zeuthen-flowed color-electric correlator on a 64 ×
16 lattice at T ≈ . T c . The correlator has been tree-level-improved by using the non-flowed leading-order continuum and lattice correlators(as in Eq. A7). Right : Wilson-flowed correlator in the same setting. The correlator has been tree-level-improved by using the
Wilson-flowed leading-order continuum and lattice correlators at the same flow time. contribution of these effects by comparing the pertur-bative continuum and lattice correlator (“tree-level im-provement” [29]). This can be achieved by defining theimproved distances τ T through G normcont τ F ( τ T ) ≡ G normlatt τ F ( τ T ) (A3)which are used to shift the correlator values according to G imp τ F ( τ T ) ≡ G latt τ F ( τ T ) . (A4)Another possibility is given by multiplying the corre-lator data at their original distances with ratio of theleading-order continuum and lattice correlator: G imp τ F ( τ T ) ≡ G normcont τ F ( τ T ) G normlatt τ F ( τ T ) G latt τ F ( τ T ) . (A5)In principle it should not matter whether one improvesthe distances or the correlator values themselves. How-ever, the interpolations of the correlator data that arenecessary for the continuum extrapolation will changeslightly. The data will not be distributed evenly if oneimproves the separations, and one obtains slightly smallerlarge separations (e.g. τ T = . → τ T ≈ .
48) which isunwished-for since most of the information about the low-frequency part of the spectral function is contained at thelargest τ T . Consequently, we choose the straightforwardstrategy of Eq. A5 to perform the tree-level improvement.In the end we want to normalize our lattice measure-ments to the leading-order continuum correlator and sowe end up with G imp τ F ( τ T ) G normcont τ F ( τ T ) = G normcont τ F ( τ T ) G normlatt τ F ( τ T ) G latt τ F ( τ T ) G normcont τ F ( τ T ) = G latt τ F ( τ T ) G normlatt τ F ( τ T ) . (A6) In principle it is necessary to take both correlators onthe right side of Eq. A6 at the same flow time τ F ; how-ever, currently we have only calculated the leading-ordercolor-electric correlator under Wilson flow (Eq. A2),which is substantially different from the more compli-cated Zeuthen flow [27] that we use for the numericalsimulations.In Figure 9 we compare the Zeuthen-flowed correla-tor improved by the non-flowed leading-order correla-tors (left panel) with the
Wilson-flowed correlator im-proved by the
Wilson-flowed leading-order correlators atthe same flow time (right panel) on a coarse lattice with N τ =
16. While the behavior under flow of these twoquantities differs greatly, their naive flow-time-to-zero ex-trapolations using a linear
Ansatz (see also section IV D)differ only by about 1%, which is shown in Figure 10.This difference is mainly caused by the discretization,and we expect it to be even smaller for finer lattices,which is why we decide to use the non-flowed tree-levelimprovement on the Zeuthen-flowed measurements: G imp τ F ( τ T ) ≡ G normcont τ F = ( τ T ) G normlatt τ F = ( τ T ) G Zeuthenlatt τ F ( τ T ) . (A7)The continuum and flow-time-to-zero extrapolations areperformed on this improved correlator normalized by theleading-order continuum correlator. Appendix B: Electric field renormalization underWilson flow
In section IV B we briefly comment on the initial risingbehavior of the color-electric correlator which is visible2 .
25 0 .
30 0 .
35 0 .
40 0 .
45 0 . τ T . . . . . . G Wilsonlatt ( τ,τ F ) /G normlatt ( τ,τ F ) (cid:12)(cid:12)(cid:12) τ F → G Zeuthenlatt ( τ,τ F ) /G normlatt τ F =0( τ ) (cid:12)(cid:12)(cid:12) τ F → N τ = 16 FIG. 10. Ratio of the Zeuthen-flowed correlator with non-flowed tree-level improvement and the Wilson-flowed correla-tor with Wilson-flowed tree-level improvement on the 64 × τ T . for the second smallest distance and beyond if the sta-tistical errors are small enough. (The rising behavior isnot visible for the higher separations as here the domi-nant effect of the flow is the improvement of the signal.)Figure 11 shows the Wilson-flowed (not Zeuthen) color-electric correlator as a function of flow time on a N τ = G QED ( τ, τ F ) = ⟨ e E i ( τ, ⃗ , τ F ) e E i ( , ⃗ , τ F )⟩ . (B1)This enables us to qualitatively analyze the renormal-ization behavior of the electric field under gradient flowwith a smaller number of diagrams compared to QCD.Nonetheless, the renormalization properties should bethe same in both theories.The leading order contribution has already been de-rived in the thermal theory (Eq. A2) and so we find G LOQED ( τ, τ F ) = e ∫ K e ik τ e − τ F ˜ K ˜ K [ k + ˜ K ] . (B2)At next-to-leading order in the lattice spacing a thereare three diagrams, which arise from the expansion of theelectric field operator (Figure 12), the action (Figure 13)and the Wilson flow equation (Figure 14).The tadpole diagram (Figure 12) is a lattice artifactand describes the mixing between the electric field op-erator E and higher-order operators such as E which τ F T . . . . . . . G latt τ ( τ F ) G normlatt τ,τ F =0 τ T =0.5000.4380.3750.312 0.2500.1880.1250.0620.0 0.05 0.08 0.10 0.12 0.13 0.14 0.15 √ τ F T FIG. 11. Non-perturbative Wilson-flowed color-electric corre-lator with Wilson-flowed tree-level improvement as a functionof flow time on a 64 ×
16 lattice at 1 . T C . are induced by the lattice (lattice tadpole corrections).The contribution of this unphysical type of diagram van-ishes in the continuum limit as there is no equivalentdiagram in the continuum. However, for non-continuum-extrapolated calculations the contribution of the diagramneeds to be considered and the interpretation is alreadywell understood and we refer to [42] for a detailed anal-ysis. One of the electric fields is expanded up to O( a ) leading to δ ( a ) G NLOQED ( τ, τ F )= − a e ⟨( G i ( X, τ F )) G i ( X ′ , τ F )⟩∣ X ′ = , ⃗ x = = − a e LO ( τ, τ F ) ∫ Q e − τ F ˜ Q ˜ Q [ q + ˜ Q ]= − a e LO ( τ, τ F ) ∫ Q e − τ F ˜ Q ·„„„„„„„„„„„„„„„„„„„„„„„‚„„„„„„„„„„„„„„„„„„„„„„„¶ ≡I ( τ F ) . (B3) kq FIG. 12. Next-to-leading order diagram of the expansion ofthe color-electric field operator. On the left the diagram isdrawn in our notation; below the horizontal dotted line arethe non-flowed propagators and vertices, and the vertical axisabove the dotted line represents evolution in flow time. Onthe right side the same diagram is drawn in the notation intro-duced in [4]. We will denote this diagram to be contribution ( a ) . kq FIG. 13. Next-to-leading order diagram of the expansion inthe action. We call this diagram ( b ) . This diagram is a pure multiplicative correction to theleading order and describes the tadpole renormalizationof the operators. The integral can be calculated analyti-cally yielding I ( τ F ) = ∫ Q e − τ F ˜ Q = e − τ F / a ( I ( τ F / a ) a ) , (B4)where I denotes the modified Bessel function of thefirst kind. The term − / a I ( τ F ) is shown in Fig-ure 15. From the figure it is clear that for zero flow timethe diagram has a significant correction to the leading-order behavior of the correlator. However, if more flowis applied, this contribution vanishes. One can inter-pret this behavior as follows: the lattice introduces high-dimension operators with coefficients containing powersof a : E → E + a E . The flow smears out the fields andintroduces a smooth momentum-cutoff, so that the ex-pectation value of such high-dimension operators is sup-pressed by an inverse power of the flow time; the effects a E are parametrically of order a / τ F2 .Another NLO effect arises from interaction terms inthe expansion of the action. The lowest order interactionis a four-point vertex. The corresponding diagram is thefirst diagram contributing to the photon self-energy up tothis order in a and is shown in Figure 13. It has the samemeaning as in the continuum: the renormalization of thecoupling. Since the loop momentum lives at zero flowtime, the corresponding correction is expected to be in-dependent of the flow time. Consequently, the correlatoras a function of gradient flow takes the form δ ( b ) G NLOQED ( τ, τ F )= − a e ∫ y ⟨ G i ( X, τ F ) F αβ ( Y ) G i ( X ′ , τ F )⟩∣ X ′ = , ⃗ x = = − a e ( τ, τ F )I , (B5)where the integral I is defined as I = ∫ Q = a − . This shift is exactly the leading-order contribution to thedifference between the bare and renormalized coupling.The last diagram (Figure 14) emerges from the non-linear part of the flow equation leading to a four-point kq k
FIG. 14. Next-to-leading order diagram of the expansion ofthe flow equation. Only tree structures are possible for τ F > ( c ) . vertex at an intermediate flow time τ ′ . It is a pure gra-dient flow correction to the photon self-energy. In QEDthis is a pure lattice artifact; in QCD similar terms alsoappear in the continuum theory but the lattice introducesa distinct gauge-invariant subclass of effects with no con-tinuum analogue. The unique NLO correction in latticeQED is δ ( c ) G NLOQED ( τ, τ F )= − a e ∫ Y τ F ∫ dτ ′ ⟨ G i ( X, τ F ) ∂ yα L β ( Y, τ ′ )× G αβ ( Y, τ ′ ) G i ( x ′ , τ F )⟩∣ X ′ = , x = = a e ∫ K e ik τ e − τ F ˜ K ˜ K ˜ K [ k + ˜ K ] ∫ Q ( e − τ F ˜ Q − ) ˜ Q ·„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„‚„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„„¶ ≡I ( τ F ) = a e ∫ K e ik τ e − τ F ˜ K [ k + ˜ K ] I ( τ F ) , (B6)where L β ( Y, τ ′ ) is the Lagrange multiplier field neededto introduce the Wilson flow equation in its Lagrangianform, see [4] for more details. We note that in this casethe integral of the leading-order contribution is not recov-ered, instead an additional power of ˜ K emerges. Sincethis diagram represents a pure Wilson flow effect, it isnot surprising that I ( τ F ) vanishes for zero flow time.The integral as a function of the flow time is shown inFigure 15. From the figure we can learn that the contri-bution of this diagram asymptotically reaches a constantvalue for large flow times. This is essentially a shift inthe flow time that one can absorb into the definition ofthe flow time. To see this, we consider the case thatmore flow loops are attached to the propagator connect-ing the electric fields. These additional loops will lead toa correction of the formLO ( τ, τ F ) ( + c ( τ F ) a ˜ K (B7) + c ( τ F ) a K + c ( τ F ) a K + . . . ) , where the coefficient functions c n ( τ F ) are determined bythe correlator at a given order in the lattice spacing a .The sum of all diagrams with iterated simple tadpoles of4 . . . . τ F /a − . − . − . − . − . − . − . − . . − a I ( τ F ) a I ( τ F ) FIG. 15. Two of the integral terms of the next-to-leading-order lattice perturbation theory calculation of the electriccorrelator in QED. The term involving I represents the tad-pole renormalization as a function of flow time and is respon-sible for the initial rising behavior of the correlator. The terminvolving I is a pure Wilson flow effect and vanishes for zeroflow time. the sort shown in Figure 15 lead to c n ( τ F ) = ( c ( τ F )) n ,and the corrections resum into the form e + ˜ K f ( τ F ) , whichcan be interpreted as a finite order- a renormalization ofthe flow time. This effect has been observed before inprevious lattice studies in [43, 44].Combining the contributions of the diagrams, we cannow write down the full correlation function at next-to-leading order: G NLOQED ( τ, τ F ) = ∫ K e ik τ e − τ F ˜ K ˜ K [ k + ˜ K ] (B8) × [ e + a e (− I − I ( τ F ) + I ( τ F ) ˜ K )] ..