Coulomb Artifacts and Breakdown of Perturbative Matching in Lattice NRQCD

By studying an explicit analytical solution of the Coulomb problem on the lattice we demonstrate a breakdown of perturbative matching for the description of the Coulomb artifacts in lattice NRQCD, which leads to a large systematic error in the predictions for the heavy quarkonium spectrum. The breakdown is a result of a fine interplay between the short and long distance effects specific to the lattice regularization of NRQCD. We show how the problem can be solved within the Schrodinger matching procedure.


Introduction
The lattice simulations within the effective theory of nonrelativistic QCD (NRQCD) [1,2,3,4] has developed into one of the most powerful tools for the theoretical analysis of heavy quarkonium properties [5]. This method is entirely based on first principles, allows for simultaneous treatment of dynamical heavy and light quarks, and gives a systematic account of the long distance nonperturbative effects of the strong interaction.
A crucial part of the approach is the matching of lattice NRQCD to the full theory of continuum QCD, which properly takes into account the effect of the hard relativistic modes. For a long time this procedure was thought to be well understood. Recently, however, a problem of the standard perturbative matching in the analysis of the heavy quark-antiquark bound states on the lattice has been pointed out [6]. The problem is related to the description of the lattice Coulomb artifacts resulting from the effect of the space discretization on the Coulomb bound state dynamics. The Coulomb artifacts appear as powers of a dimensionless combination α s am q of the strong coupling constant α s , lattice spacing a, and the heavy quark mass m q in the parameters of the bound states evaluated on the lattice. This dependence on the lattice spacing should be cancelled in the final result for the physical quarkonium spectrum through the matching procedure and an inconsistent treatment of the artifacts may lead to a large systematic error of the lattice NRQCD predictions. In the present paper, by studying an explicit analytical solution of the Schrödinger equation on the lattice, we show that the standard finite-order perturbative matching is insufficient in dealing with the Coulomb artifacts, confirming the result of numerical analysis [6]. This appears counterintuitive since the lattice regularization is usually associated with a momentum cutoff at the scale 1/a much larger than the scale α s m q of Coulomb dynamics and the corresponding short-distance effects are supposed to be systematically described by the standard matching procedure. We find that the failure of perturbative matching is a consequence of a fine interplay between the short and long distance effects specific to the lattice regularization of NRQCD. The paper is organized as follows. In the next section we consider the case of the four-quark operators in the NRQCD Lagrangian where the perturbative matching breaks down in one loop, and show how the problem can be related to the properties of the Coulomb wave function on the lattice within the method of "Schrödinger matching". In Sects. 3 and 4 we describe an exact analytical solution of the Schrödinger equation with Coulomb potential on a spatial lattice, which is used to model the Coulomb artifacts of the real lattice NRQCD simulations without the expansion in α s , and explain the mechanism of the perturbative matching breakdown. In Sect. 5 we discuss the application of our result to the analysis of the bottomonium spectrum in lattice NRQCD.

Perturbative and Schrödinger matching of four-quark operators
The effect under consideration is characteristic to the four-quark operators contribution to the NRQCD Lagrangian where ψ (χ c ) are the nonrelativistic Pauli spinors of quark (antiquark) field, Γ i is a matrix in color and spinor space and C i is the corresponding Wilson coefficient. The coefficients C i vanish in the Born approximation and is determined order by order in α s by equating the one-particle irreducible quark-antiquark scattering amplitudes in QCD and NRQCD in a given order in heavy quark velocity v. The Coulomb artifacts are associated with the static Coulomb interaction between the nonrelativistic quark and antiquark. The nonrelativistic Coulomb dynamics is not sensitive to the high momentum region and, to the leading order in v, the result for such a contribution obtained within continuum QCD and NRQCD coincide.
Thus the matching is determined by the difference between continuum and lattice NRQCD amplitudes. In the one-loop approximation the relevant NRQCD diagram is given in Fig. 1, diagram as well as for the diagrams with any number of the Coulomb exchanges can be incorporated into the "Schrödinger matching" for the value of the wave function at the origin introduced in Ref. [6]. The corresponding correction to the energy level is proportional to where ψ c (ψ l ) stands for the Coulomb wave function computed in the continuum (on the lattice). One may suggest that for α s m q ≪ 1/a the two procedures should give equivalent results and the Coulomb artifact contribution to C i should reproduce Eq. (3) in a given order of the expansion in α s . However in Ref. [6] it has been demonstrated by numerical analysis becomes with the boundary condition R(0) = 0. Throughout this section we use the Coulomb units, with the reduced Bohr radius r B = (C F α s m q /2) −1 being the unit of length. The analytical solution of the eigenstate problem Eq. (4) is known and the corresponding eigenfunctions for discrete and continuum spectrum are found in terms of a hypergeometric function [7]. For the ground state the solution takes a particulary simple form where E c = − 1 2 and ψ c (r) = 1 √ π e −r are the well known continuum results for the ground state energy and wave function, and the normalization condition reads The wave function in the momentum space is defined by the Fourier transform and readsψ whereψ c (p) = 8 √ π (p 2 +1) 2 is the continuum result in Coulomb units. A peculiar feature of the Fourier transform is that ψ l (0) does not contribute to Eq. (8) due to the n 2 factor. As a consequence the inverse Fourier transform ψ l (n) = 1 2π 2 π a 0 dp p 2ψ l (p) give the correct value of ψ l (n) only for the sites with n = 0. For n = 0 the corresponding takes the value while the original solution in the coordinate space is in full agreement with the results of the numerical analysis [6]. Though which in the limit a → 0 recovers the property of the continuum solution At the same time if Eq. (12) is used to define the value of the wave function at the origin we get which violates the continuum result at O(1).
Finally we should note that while the linear term in Eq.

Breakdown of perturbative matching
Let us now consider the Coulomb contribution of the diagram Fig. 1 to a Wilson coefficient C i . As it has been pointed out this contribution is given by the difference of the continuum and lattice NRQCD one-loop expressions and can be computed in the kinematics where the heavy quarks are at rest and on the threshold. In the given order of the nonrelativistic expansion it reads where the integration over the time component of the virtual momentum has been performed by taking the residue of the heavy quark propagator, an auxiliary infrared cutoff λ is introduced to regularize the integrals at small momentum, and we switch back to the standard units. The continuum Coulomb gluon and heavy quark propagators correspond to the leading order NRQCD action. Their lattice counterparts D l (p) = sin(ap) ap (a/2) 2 sin 2 (ap/2) , G l (p) = m q (a/2) 2 sin 2 (ap/2) exactly correspond to the discretization of the Schrödinger equation (4). After integration we get Keeping the most singular terms in the limit r → 0 we get which reproduces Eq. (15) in Coulomb units. The first (second) term in Eq. (22) corresponds to the kinetic (potential) energy contribution. Hence at the origin these terms should be considered on an equal footing while the standard matching treats the Coulomb potential as a perturbation. Evidently the above mechanism of the perturbative matching breakdown is specific to the lattice regularization and the contact interaction in the NRQCD Lagrangian.
We can regularize the contact interaction by separating the quark and antiquark fields with a small spatial interval r 0 ∼ a and take the limit r 0 → 0 after the matching is done. This introduces an additional sin(r 0 p)/(r 0 p) factor into the integrands of Eq. (17). The lattice integral then reads m q π a λ dp p 2 sin(r 0 p) r 0 p sin(ap) ap (a/2) sin(ap/2) As we have shown above, the perturbative matching of the four-quark operators does not correctly account for the Coulomb artifacts and, for a ∼ 1/(vm b ) and v ∼ α s , results in O(1) error in the prediction for the spectrum. For example, the one-loop matching of the spin-flip four-quark operator with the spurious linear Coulomb artifact gives the value of the bottomonium hyperfine splitting [12], which overshoots the predictions of perturbative QCD [13] by almost a factor of two, in clear conflict with the general understanding of the heavy quarkonium dynamics.
In practice the effect of the lattice artifacts is reduced by numerical extrapolation of the data to a = 0 [5,11]. The extrapolation below a ∼ 1/m b in this case is justified because for the typical values of lattice spacing the numerical effect of the 1/(am b ) n terms on the data points is small. This extrapolation effectively removes all the lattice artifacts including the (aΛ QCD ) n terms associated with the effect of lattice regularization on the dynamics at the confinement scale Λ QCD . 4 The problem of the perturbative matching breakdown, however, is not fully fixed by this procedure. Since the radiatively improved lattice result is supposed to be free of linear artifacts, the extrapolation is performed through a constrained fit of the data points by a polynomial in a with vanishing linear term. Since the bare lattice data are free of the linear Coulomb artifacts the one-loop perturbative matching in fact introduces a linear dependence of the radiatively improved result on a, which leads to a systematic error of the fit. For example, in the analysis of the bottomonium hyperfine splitting [11] based on the one-loop perturbative matching [12,14] this error exceeds 10% of the total result and is beyond the estimated uncertainty interval of the radiatively improved O(v 6 ) lattice NRQCD analysis [6]. At the same time the perturbative matching can be used for the self-consistent analysis of the quarkonium spectrum within the above extrapolation scheme if the Coulomb artifacts are removed from the Wilson coefficients by means of the asymptotic expansion [15] or a numerical fit [6].