{'fraction_non_alphanumeric': 0.04720305121685434, 'fraction_numerical': 0.015655648383581546, 'mean_word_length': 4.043601722084822, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 20, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
New lattice action for heavy quarks

We extend the Fermilab method for heavy quarks to include interactions of dimensions 6 and 7 in the action. There are, in general, many new interactions, but we carry out the calculations needed to match the lattice action to continuum QCD at the tree level, finding six non-zero couplings. Using the heavy-quark theory of cutoff effects, we estimate how large the remaining discretization errors are. We find that our tree-level matching, augmented with one-loop matching of the dimension-five interactions, can bring these errors below 1%, at currently available lattice spacings.


I. INTRODUCTION
An important application of lattice gauge theory is to calculate hadronic matrix elements relevant to experiments in flavor physics. With recent advances in lattice calculations with n f = 2 + 1 flavors of dynamical quarks [1][2][3][4], we now have an exciting prospect of genuine QCD calculations. To match the experimental uncertainty, available now or in the short term, it is essential to control all other sources of theoretical uncertainty as well as possible. An attractive target is to reduce the uncertainty, from any given source, to 1-2%. This target will be hard to hit if one relies on increases in computer power alone: methodological improvements are needed too.
Many of the important processes are electroweak transitions of heavy charmed or b-flavored quarks. A particular challenge stems from heavy-quark discretization effects, because m Q a ≪ 1.
The key to meeting the challenge is to observe that heavy quarks are non-relativistic in the rest frame of the containing hadron [5,6]. The scale of the heavy-quark mass, m Q , can (and should) be separated from the soft scales inside the hadron and treated with an effective field theory instead of computer simulation. Even so, at available lattice spacings [1], many calculations of D-meson (B-meson) properties suffer from a discretization error of around 7% (5%) [2,3]. Thus, it makes sense to develop a more accurate discretization.
In this paper we extend the accuracy of the "Fermilab" method for heavy quarks [7] to include in the lattice action all interactions of dimension six. We also include certain interactions of dimension seven. Because heavy quarks are non-relativistic, they are commensurate with related dimension-6 terms, in the power counting of heavy-quark effective theory (HQET) for heavy-light hadrons [5] or non-relativistic QCD (NRQCD) for quarkonium [6].
The Fermilab method starts with Wilson fermions [8] and the clover action [9]. With these actions lattice spacing effects are bounded for large m Q a, thanks to heavy-quark symmetry. They can be reduced systematically by allowing an asymmetry between spatial and temporal interactions. Asymmetry in the lattice action compensates for the non-relativistic kinematics, enabling a relativistic description through the Symanzik effective field theory [10]. Alternatively, one may interpret Wilson fermions non-relativistically from the outset [7], and set up the improvement program matching lattice gauge theory and continuum QCD to each other through HQET and NRQCD [11,12]. The Symanzik description makes it possible to design a lattice action that behaves smoothly as m Q a → 0, converging to the universal continuum limit. The HQET description, on the other hand, makes semiquantitative estimates of discretization errors more transparent.
The new action introduced below has nineteen bilinear interactions beyond those of the asymmetric version of the clover action, as well as many four-quark interactions. Several of these couplings are redundant, and many more vanish when matching to continuum QCD at the tree level. We study semiquantitatively how many of the new operators are needed to achieve 1-2% accuracy. We find, in the end, that only six new interactions are essential for such accuracy. The action is designed with some flexibility, so that one may choose the computationally least costly version of the action. This paper is organized as follows. Section II considers the description of lattice gauge theory via continuum effective field theories. Then, in some detail, we identify a full set of operators describing heavy-quark discretization effects. We then determine how many of these are redundant, and which redundant directions should be used to preserve the good high-mass behavior. We have two goals in this analysis. One is to design the new, more highly improved, action; for this step a Symanzik-like description is more helpful, and the resulting action is given in Sec. III. The other is to estimate the discretization errors of the new action; here the HQET and NRQCD descriptions are more useful. To make error estimates, and to use the new action in numerical work, we need matching calculations; they are in Sec. IV. Our error estimates are in Sec. V. Section VI concludes. Some of the material is technical and appears in appendices: Feynman rules needed for the matching calculation are in Appendix A; some details of the Compton scattering amplitude used for matching are in Appendix B; a discussion of improvement of the gauge action on anisotropic lattices (which one needs only if the heavy quarks are not quenched) is in Appendix C. Some of these results have been reported earlier [13].

II. EFFECTIVE FIELD THEORY
In this section we discuss how to understand and control discretization effects using effective field theories. We start with a brief overview, focusing on issues that arise for heavy quarks, those with mass m Q ≫ Λ. For more details, the reader may consult earlier work [7,11,12,14,15] or a pedagogical review [16]. Here we catalog all interactions of dimension 6 and also certain interactions of dimension 7 that, for heavy quarks, are of comparable size when m Q a ≪ 1.

A. Overview
Cutoff effects in lattice field theories are most elegantly studied with continuum effective field theories. The idea originated with Symanzik [10] and was extended to gluons and light quarks by Weisz and collaborators [9,[17][18][19]. One develops a relationship where .
= means that the two Lagrangians generate the same on-shell spectrum and matrix elements. The lattice itself regulates the ultraviolet behavior of the underlying (lattice) theory L lat . On the other hand, a continuum scheme, which does not need to be specified in detail, regulates (and renormalizes) the ultraviolet behavior of the effective theory L Sym .
In lattice QCD (with Wilson fermions), the local effective Lagrangian (LEL) is where g 2 and m f are the gauge coupling and quark mass (of flavor f ), renormalized at scale µ a −1 . The (continuum) QCD Lagrangian appears as the first two terms. The sum consists of higher dimension operators L i , multiplied by short-distance coefficients K i . These terms describe cutoff effects. The short-distance coefficients depend on the renormalization point and on the couplings, including couplings c j of improvement terms in L lat . Equation (2.2) is fairly wellestablished to all orders in perturbation theory [20,21] and believed to hold non-perturbatively as well. If a is small enough, the terms L i may be treated as operator insertions, leading to a description of lattice gauge theory as "QCD + small corrections". In heavy-quark physics m Q ≫ Λ, where Λ is the QCD scale, so one is led to consider what happens when m Q a ≪ 1. The short-distance coefficients depend explicitly on the mass. Time derivatives of heavy-quark or heavy-antiquark fields in the L i also generate mass dependence of observables. With field redefinitions-or, equivalently, with the equations of motion-these time derivatives can be eliminated. Focusing on a single heavy flavor Q, the result of these manipulations is [7,14,15] L Sym = · · · −Q γ 4 D 4 + m 1 + m 1 m 2 γ · D Q + i a dimL i −4K i (g 2 , m 2 a; µa)L i , (2.3) where the ellipsis denotes the unaltered LEL for gluons and light quarks. By construction theL i do not have any time derivatives acting on quarks or antiquarks. The advantage of Eq. (2.3) is that all dependence on the heavy-quark mass is in the shortdistance coefficients m 1 , m 1 /m 2 , andK i (m 2 a). Matrix elements of theL i generate soft scales. The heavy-quark symmetry of Wilson quarks (with either the Wilson [8] or Sheikholeslami-Wohlert [9] actions) guarantees that the coefficientsK i (m 2 a) are bounded for all m 2 a. This feature can be preserved by improving the lattice Lagrangian with discretizations of theL i , thereby avoiding higher time derivatives [7,11]. For such improved actions, Eq. (2.3) neatly isolates the potentially most serious problem of heavy quarks into the deviation of the coefficient m 1 /m 2 from 1.
Fortunately, the problem can be circumvented in two simple ways. One is a Wilson-like action with two hopping parameters [7], tuned so that m 1 = m 2 . Then Eq. (2.3) once again takes the form "QCD + small corrections". The new lattice action introduced in Sec. III has two hopping parameters for this reason.
Another solution is to interpret Wilson fermions in a non-relativistic framework. One can replace the Symanzik description with one using a non-relativistic effective field theory for the quarks (and antiquarks) [11]. For the leadingQ-Q term in Eq. (2.3) (2.4) where z B is a matching coefficient, and h (+) is a heavy-quark field satisfying h (+) = +γ 4 h (+) . Another set of terms appears for the antiquark, with field h (−) satisfying h (−) = −γ 4 h (−) . The non-relativistic effective theory conserves heavy quarks and heavy antiquarks separately. As a consequence, the rest mass m 1 has no effect on mass splittings and matrix elements. 1 For lattice gauge theory this implies that the bare quark mass (or hopping parameter) should not be adjusted via m 1 . Instead, the bare mass should be adjusted to normalize the kinetic energy D 2 /2m 2 .
One can develop the non-relativistic effective theory for the lattice artifactsL i by using heavyquark fields instead of Dirac quark fields [11]. Higher-dimension operators in the heavy-quark theory receive contributions from the expansions of Eq. (2.4) and of theL i . Coalescing the coefficients of like operators obtains a description of lattice gauge theory with heavy quarks where the operators O i on the right-hand side are those of a (continuum) heavy-quark effective theory, of dimension 5 and higher, built out of heavy-quark fields h (±) , gluons, and light quarks.
(The leading ellipsis denotes term for the gluons and light quarks only.) The C i are short-distance coefficients, which depend on g 2 , the heavy-quark mass, the ratio of short distances m 2 a, and also all couplings c j in the lattice action. The logic and structure is the same as the non-relativistic description of QCD, Thus, improvement of lattice gauge theory is attained by adjusting couplings c j until C lat i (c j )−C cont i vanishes (identically, or perhaps to some accuracy) for the first several O i . 1 A simple proof can be found in Ref. [11].
It does not matter whether one carries out the improvement program by adjustingK i (c j ) = 0 or C lat i (c j ) = C cont i [12]. The results for the c j are the same, provided one identifies m Q with m 2 . The matching assumes that pa ≪ 1, but at the same time m 2 a ≪ 1. One is thus led to non-relativistic kinematics (p/m 2 ≪ 1) in the matching calculation, where both descriptions-Eqs. (2.3) and (2.5)-are valid. Kinematics are encoded into the operatorsL i or O i and are not transferred to the short-distance coefficients. Hence, kinematics cannot influence matching conditions on the c j . In particular, when indeed m 2 a ≪ 1 (which may be impractical, but is conceivable theoretically) relativistic kinematics (p ∼ m 2 ) are possible, and it follows from the Symanzik effective field theory that the solution ofK i (c j ) = 0 yields the same c j for both relativistic and non-relativistic kinematics.

B. Quark bilinears in the LEL
In the rest of this section we construct the LEL appropriate to heavy quarks. The two main steps are first to list all of the L i that can appear, and second to decide which should be considered redundant. In part it is a generalization of the dimension-6 analysis of Ref. [9] to the case without axis-interchange symmetry. At dimension 6 there are quark bilinears, four-quark interactions, and interactions that contain only the gauge field. We shall start with the bilinears and turn to the others further below. In each case, we first consider complete lists of operators, and then consider which can be chosen to be redundant. Table I contains a list of all quark bilinears through dimension 6 that can appear in the effective Lagrangian. The second column contains interactions that respect axis-interchange symmetry; the fourth column contains the extension to the case without axis-interchange symmetry. The meaning of the other columns is explained below. Covariant derivatives act on all fields to the right, (2.7) This notation is convenient for the interactions with commutators and anti-commutators. To arrive at the lists we exploit identities such as (2.10) Some interactions are omitted, because the underlying lattice gauge theory is invariant under cubic rotations, spatial inversion, time reflection, and charge conjugation. 2 The fourth column is arranged so that its entries are part of the corresponding interactions in the second column. It is easy to show that the list is complete, by writing out all independent ways to have three covariant derivatives, expressing the E and B fields as anti-commutators of covariant derivatives. One finds 11 possibilities, and then one can use identities to manipulate this list to that given in the fourth column of Table I.
The LEL contains several redundant directions. The equation of motion of the leading LEL plays a key role in specifying which operator insertions may be considered redundant. Let us assume, for the moment, that m 1 = m 2 , so that the equation of motion in the Symanzik LEL is the Dirac equation. Below we shall use the non-relativistic effective field theory to address the case The quark fields are integration variables in a functional integral, so an equally valid description is obtained by changing variables (2.13) and similarly forJ with separate parametersε i ,δ i , andθ i . If the δ parameters (and ϑ 2 ,θ 2 ) vanish, then J andJ preserve invariance under interchange of all four axes. One can propagate the change of variables to the LEL, and trace which coefficients of dimensions 5 and 6 are shifted by amounts proportional to the parameters in J andJ. To avoid generating terms that violate charge conjugation one choosesε i = +ε i ,δ i = +δ i ,θ 2 = −ϑ 2 . We then see that there are two redundant directions at dimension 5, and five at dimension 6. That means that two couplings in the dimension-5 lattice action may be set by convenience, and five in the dimension-6 lattice action. The third and fifth columns show the correspondence between parameters in the change of variables and the interactions that we choose to be redundant. As expected from general arguments [7,14,15], all interactions in which γ 4 D 4 acts on Q or (after integration by parts)Q are redundant.
There is quite a bit of freedom here. One could choose ε F to eliminateQ[D 4 , γ · E]Q = Q{γ 4 D 4 , α · E}Q instead ofQ{γ · D, α · E}Q. But the former is suppressed, relative to the latter, in heavy-quark systems. Moreover, in HQET and NRQCD one has which mean thatQα · EQ andQ{γ · D, α · E}Q generate nearly the same effects in heavy-quark systems. Thus, we prefer to takeQ{γ · D, α · E}Q to be redundant.
To understand the general pattern of redundant interactions, let us introduce some notation. Let B (E) be a combination of gauge fields, derivatives, and Dirac matrices that commutes (anticommutes) with γ 4 . An example of B (E) is iΣ · B (α · E). Also, let us write B ± (and E ± ) when QB ± Q (orQE ± Q) has charge conjugation ±1.

C. Power counting
The small corrections of an effective field theory are small, because the product of the shortdistance coefficients and the operators yield a ratio of a short-distance scale to a long-distance scale. For light quarks in the Symanzik effective field theory, the essential ratio is a/Λ −1 = Λa, and dimensional analysis reveals the power of Λa to which any contribution is suppressed. In particular, Band E-type interactions of the same dimension are equally important.
For heavy quarks the physics is different, because m −1 Q is a short distance. The ratio a/m −1 Q = m Q a should not be taken commensurate with Λa [7]. Instead, interactions should be classified in a way that brings out the physics. It is natural to turn to HQET and NRQCD. Let us start with heavy-light hadrons and HQET. E-type interactions of given dimension are Λ/m Q times smaller than B-type interactions of the same dimension. Because Λ/m Q ≪ 1 and Λa ≪ 1, it makes sense to count powers of λ, where λ is either of the small parameters [11,12,15] λ ∼ aΛ, Λ/m Q . Here n Γ = 0 or 1 for interactions of the formQB + Q orQE + Q, respectively. The sixth column of Table I (labelled HQET) shows the suppression of each interaction, relative to the (leading) contribution from the light degrees of freedom. In the following we call the power counting for heavy-light hadrons, based on Eq. (2.17), "HQET power counting." Now let us recall how to classify interactions in quarkonium according to the power of the relative internal velocity, υ. Because color source and sink are both non-relativistic, chromoelectric fields carry a power of υ 3 , and chromomagnetic fields a power of υ 4 [22]. E-type interactions are suppressed by a power of p/m Q = υ, analogously to their suppression in heavy-light hadrons. Thus, bilinears are suppressed by υ t , where now t = dim L − 3 + n E + 2n B + n Γ , (2.18) and n E (n B ) is the number of chromoelectric (chromomagnetic) fields. The seventh column of Table I (labelled NRQCD) shows the suppression of each interaction. In the following we call the power counting for quarkonium, based on Eq. (2.18), "NRQCD power counting." Glancing down the sixth and seventh column of Table I, one sees several terms of order λ 3 and υ 6 , from Eqs. (2.17) and (2.18) one realizes that some dimension-7 interactions are of the same order. They are listed in Table II. There are two interactions with four derivatives, six with the chromomagnetic field and two derivatives, and four with two E or two B fields. A third combination of four derivatives is omitted, using the identity Other dimension-7 operators carry power λ 4 in HQET power counting, or υ 8 (or higher) in NRQCD power counting. Five combinations are redundant (as shown), and we shall see below how they and the others arise in matching calculations.
The (d, n Γ ) = (7, 1) operatorQ{D 2 , α · E}Q and several (d, n Γ ) = (8, 0) operators, all with n E = 1 and n D + n Γ = 3, have NRQCD power-counting υ 6 . Reference [22] includes spindependent ones, to obtain the next-to-leading corrections to spin-dependent mass splittings. We have not included these operators in our analysis, but a straightforward extension of the matching calculation in Sec. IV B 1 would suffice to determine their couplings.
Although this description of cutoff effects is somewhat cumbersome, it provides a valuable foundation for our new action, given in Sec. III. To obtain the new action, we simply discretize the interactions in Tables I and II, except those with higher time derivatives. The discretization of Qγ · DQ is needed to obtain a lattice action that behaves smoothly as m Q a → 0 [7], reproducing the universal continuum limit of QCD. Similarly, discretizations of the E-type interactions, such asQα · EQ andQ{γ · D, D 2 }Q, are needed to retain that feature here.

Dim
Without axis-interchange symmetry

D. Heavy-quark description
For understanding the size of heavy-quark discretization effects, it is simpler to switch to a non-relativistic description. (When m 1 = m 2 , it is also necessary to see the connection to QCD.) The list of interactions is much shorter, because the constraint γ 4 h (±) = ±h (±) removes the Etype interactions. It is given in Table III, including the dimension-7 interactions related to those in Table II. Also, fewer changes of the field variables are possible: where now and similarly forJ. To avoid C-odd interactions, one should choose equal parameters in J andJ . Thus, there are four redundant directions of interest-all with time derivatives of the (anti-)quark field. In the end, just as many non-redundant interactions remain as in the Symanzik description. The heavy-quark description provides a good way to estimate the size of remaining discretization effects, as in Sec. V.

E. Gauge-field and four-quark interactions in the LEL
We now turn to interactions in the gauge sector of the LEL, and also to four-quark interactions. The two are connected when one considers on-shell improvement, because in quark-quark scattering short-distance gluon exchange generates the same behavior as four-quark contact interactions. Here we give a cursory sketch of the gauge action. Then we consider the four-quark interactions, including details mostly for completeness. In practice (see Sec. V), we find the four-quark corrections to be smaller than those of the bilinear interactions analyzed in the preceding subsection.
The gauge sector of the LEL is the same as for anisotropic lattices, where one adjusts the action so that the temporal lattice spacing a t differs from the spatial lattice spacing a s . The short-distance coefficients are different; here asymmetry between spatial and temporal gauge couplings arise only from heavy-quark loops. Improved anisotropic actions have been discussed in the literature [23], but full details remain unpublished [24]. We present the details in Appendix C.
We are most concerned here with effects that survive on shell, so we study here the possible changes of variables for the gauge field. With axis-interchange symmetry one has [9,19] with a color-adjoint vector-current term for each flavor f of quark (heavy or light). The appearance of g 2 multiplying the currents is a convenient normalization convention. When one now considers giving up axis-interchange symmetry, one has which reduce to Eq. (2.22) when the δs vanish. For a moment, let us set ε Jf = δ Jf = 0 in Eqs. (2.23) and (2.24), and focus on the gauge fields alone. As discussed in Appendix C, there are eight independent gauge-field interactions that arise at dimension six. There are three independent ways-parametrized by ε A , δ A , and δ E -to transform the gauge field, yielding three redundant directions. Similarly, there are eight distinct classes of six-link loops, shown in Fig. 1, that can be used in an improved lattice gauge action. In Appendix C, we show that three of them-all three classes of "bent rectangles" in the bottom row of Fig. 1-may be omitted from an on-shell improved gauge action.
The transformations involving the currentsq f γ µ t a q f are more interesting. They shift the LEL [cf. Eq. (2.2)] by where the derivatives act only on the gauge fields. The size of these shifts-of order g 2 for fourquark operators and of order g 0 for bilinears-is commensurate with the respective terms that already appear in L Sym . Thus, the 2n f parameters ε Jf and δ Jf could be used to eliminate bilinears or four-quark operators. For simulations it is preferable to remove the latter, We now list the dimension-six four-quark interactions in the LEL. For a single flavor, the complete list is in Table IV, which also indicates that the current-current interactions are redundant. Interactions with the color structure (qΓq) 2 may be omitted, because they can be related to those listed through Fierz rearrangement of the fields.
When considering several flavors of quark, we must keep track of flavor indices as well as color and Dirac indices. The Fierz problem becomes more intricate, and we shall find that color-singlet and color-octet structures should be maintained. Let us start with Fierz rearrangement of the Dirac indices. The four-quark terms in the LEL take the form where K X denotes short-distance coefficients, the Greek (Latin) indices label color (flavor), F is (qt a q) 2 (Qt a Q) 2 (qγ 5 t a q) 2 ( the Fierz rearrangement matrix (with F 2 = 1), and the minus sign comes from anti-commutation of the fermion fields. Equation (2.26) leaves the flavor and color indices uncontracted, but to get terms in the LEL, the color indices must be contracted (one way or another), and the flavor labels must yield a flavor-neutral interaction. Without loss, we can choose the side of Eq. (2.26) such that the Dirac matrices contract quark fields of the same flavor. Then one can use Fierz identities so that the color indices are contracted across the same fields as the Dirac and flavor indices. After using Fierz rearrangement to bring quarks of the same flavor next to each other, one is left with the interactions in Table V. To be concrete, we consider n l flavors of light quarks (with m q Λ) and two flavors of heavy quarks (charm and bottom). We neglect the dependence of the coefficients on the light quark masses, because four-quark interactions are already small corrections (of dimension six). In that case, the four-quark interactions can be arranged so that only the SU(n l ) flavor singlets fq f Γ X t a q f and fq f Γ X q f appear.
The parameters ε Jf and δ Jf may be used to eliminate color-octet current-current interactions. For each heavy flavor, one finds (Qγ 4 t a Q) 2 and i (Qγ i t a Q) 2 to be redundant. For light quarks, we may neglect the differences in the mass, so they have common parameters, and the flavorsinglet combination ( fq f γ µ t a q f ) 2 is redundant. For the light flavors, our list of operators is a Fierz rearrangement of the list in Ref. [9].
The leading HQET power counting for heavy-light four-quark operators follows from dimensional analysis and Eq. (2.17): λ 2+n Γ , just as if the light-quark part were replaced by three deriva-tives. Heavy-heavy four-quark operators will be suppressed, once matrix elements are taken, by a heavy-quark loop, leading to g 2 λ 4+n Γ .
In quarkonium, the size of heavy-light four-quark operators follows similarly from Eq. (2.18): υ 3+n Γ . The valence heavy-heavy operators are more interesting. They must contain two contributions, one to improve t-channel gluon exchange, and another to improve s-channel annihilation. The former have NRQCD power counting g 2 υ 3+n Γ ∼ υ 4+n Γ (since g 2 ∼ υ [22]). The latter are υ 2 times smaller, because the s-channel gluon is far off shell, but the Dirac-matrix suppression is now υ 1−n Γ , leading to g 2 υ 6−n Γ ∼ υ 7−n Γ in all. In practice, the s-channel contributions are suppressed further, when treated as an insertion in a color-singlet quarkonium state. At the tree level, the only color structure that can arise is the color-octet. Its matrix elements vanish in theQQ-colorsinglet Fock state of quarkonium, leaving the υ 3 -suppressedQQA color octet [25]. Color-singlet four-quark operators arise at one loop, with an additional factor of g 2 ∼ υ.

III. NEW LATTICE ACTION
In this section we introduce a new, improved lattice action for heavy quarks, designed to yield smaller discretization errors than the action in Ref. [7]. Our design is based on several lessons from the preceding section and Refs. [7,11,12]. First, it is important to preserve the natural heavy-quark symmetry of Wilson fermions, so that the coefficientsK i stay bounded for all m Q a. (This feature is spoiled in the standard improvement program designed for light quarks, which introduces several new terms that grow with m Q .) Second, the new lattice action is flexible enough to match cleanly onto both the Symanzik description and the non-relativistic description.
Let us write the action as follows where S D 2 F 2 is the improved gauge action [Eq. (C7)], S 0 is the basic Fermilab action, the S (d,n Γ ) consist of the bilinear terms added to improve the quark sector, and Sq qqq denotes four-quark interactions. S (d,n Γ ) consists of (discretizations of) interactions of dimension d, with n Γ as in the discussion of power counting, Eqs. (2.16)-(2.18). Including the interactions in S (d,1) couples "upper" and "lower" components, but allows a smooth limit a → 0. 3 Our aim is to improve the action to include all interactions of dimension six. Then the power counting requires us to include S (7,0) as well. Finally, Sq qqq consists of discretizations of four-quark operators, at dimension six, those of Table V.
The basic Fermilab action [7] is a generalization of the Wilson action [8]: We denote lattice fermions fields with ψ to distinguish them from the continuum quark fields in Sec. II. The dimension-five Wilson terms are included in S 0 to remove doubler states. The remaining dimension-five interactions are [7,9] where the notation S B and S E is from Ref. [7], and the discretizations D µ lat , △ µ lat , △ lat , B lat , E lat are defined below.
The new interactions in Eq. (3.1) introduced in this paper are All couplings in Eqs. (3.2)-(3.7) are real; explicit factors of i are fixed by reflection positivity [26] of the continuum action. Some of the improvement terms extend over more than one timeslice, so there are small violations of reflection positivity for the lattice action. We expect that the associated problems are not severe, as with the improved gauge action [27].
Equations (3.5)-(3.7) contain 19 new couplings. The convention for couplings c i , r i and z i is as follows. In matching calculations we find that couplings z i vanish at the tree level, while the couplings c i do not. Couplings r i are redundant and, for this reason, could be omitted. The analysis in Sect. II gives the number of redundant interactions, rather than the specific choices of interactions themselves. The possibilities for the dimension-7 redundant directions are as follows. One of (c 4 , c 5 , r 5 ) is redundant; we choose r 5 . Furthermore, one of (z 6 , z 7 , r 7 , r BB ), another of (z 7 , r 7 , r BB ), and another of (z 7 , r 7 , r ′ 7 , r BB ) are redundant; we choose r 7 , r ′ 7 , and r BB . But because pragmatic considerations could motivate other choices, we keep all of them in our analysis. This strategy also provides a good way for the matching calculations to verify the formal analysis of the LEL. In future numerical work, we recommend choosing r s , as usual, to solve the doubling problem (in practice r s ≥ 1). The others may be chosen to save computer time, which presumably means choosing the couplings of computationally demanding interactions to vanish.
The difference operators and fields with the subscript "lat" are taken to be where the covariant translation operators T ±ρ translate all fields to the right one site in the ±ρ direction, and multiply by the appropriate link matrix [28]. These discretizations are conventional for S 0 + S B + S E . For the new interactions, we have re-used the same ingredients. For the interactions with couplings r 5 and z ′ 7 one can consider In tree-level matching calculation, both lead to the same dependence on r 5 and z ′ 7 . Equation (3.11) has the advantage that is re-uses elements that are already defined (in a computer program, say) for the dimension-4 and -5 action. Equation (3.12) is more local, however, and may have other advantages. A FermiQCD [29] computer code of the new action indicates that Eq. (3.11) is faster [30]. This code also indicates that it is advantageous to choose the redundant directions so that one may set r 5 = r 7 = 0.
The improved gluon action S D 2 F 2 is defined in Appendix C. The four-quark action Sq qqq contains the obvious discretization of the (continuum) operators explained in Sec. II E and listed in Tables IV and V: simply substitute lattice fermion fields for the continuum fields, and assign each a real coupling. When matching to continuum QCD, the couplings in Sq qqq start at order g 2 , making them commensurate with order-g 2 matching effects in S (6,1) + S (7,0) , such as tree-level quark-quark scattering. To incorporate the four-quark action in a Monte Carlo simulation, one would introduce auxiliary fields to recover a bilinear action. In the next section we show, however, that these operators are not necessary for the target accuracy of 1-2%, so this cumbersome set-up can be avoided for now.

IV. MATCHING CONDITIONS
In this section we derive improvement conditions on the new couplings at the tree level. We calculate on-shell observables for small pa without any assumption on m Q a. We look at the energy as a function of 3-momentum, which is sensitive to c 1 , c 2 , c 4 , and z 6 . We then look at the interaction of a quark with classical background chromoelectric and chromomagnetic fields. The former is sensitive to c E , r E , and z E ; the latter to all but c EE , r EE , z EE , r BB , and z BB . To ensure that these results are compatible with the improved gauge action, we next compute the amplitude for quark-quark scattering. This step also matches the four-quark interactions, which are not written out explicitly in Sec. III. Finally, we compute the amplitude for Compton scattering to match c EE , r EE , z EE , r BB , and z BB .

A. Energy
The energy of a heavy quark on the lattice is defined through the exponential fall-off in time of the propagator. For small momentum p the energy can be written which (at fixed m 0 a) prescribes a line in the (c 1 , c 4 ) plane. As m 0 a → 0, both lines become vertical: the coefficients c 1 and c 2 of dimension-6 operators are fixed, whereas the coefficients of c 4 and z 6 dimension-7 operators are undetermined. At this stage it is tempting to choose c 4 and z 6 to be two of the redundant couplings, but below we shall see that there are better choices.

B. Background Field
To compute the interaction of a lattice quark with a continuum background field, we have to compute vertex diagrams with one gluon attached to the quark line. The Feynman rules are given in Eqs. (A23) and (A24). Our Feynman rules introduce a gauge potential via where e µ is a unit vector in the µ direction, and take the Fourier transform of the gauge field to be A background field would, however, lead to parallel transporters To use the interaction with a background classical field as a matching condition, we must compute the current J µ that couples to the background field A µ in Eq. (4.10). Current conservation requires k · J(k) = 0, (4.11) where k is the external gluon's momentum. The usual convention for A µ (k), from Eqs. (4.8) and (4.9), yields a currentĴ µ satisfyingk wherek µ = (2/a) sin(k µ a/2). One sees, therefore, that a classical gluon line with Lorentz index µ must be multiplied by One should think of n µ (k) as a wave-function factor for the external line. Its appearance has been noted previously by Weisz [17].
In the rest of this subsection we match the vertex function in lattice gauge theory with our new action to that in the continuum gauge theory. The incoming quark's momentum is p, the outgoing p ′ , and the gluon's K = p ′ − p. The current is given by (no implied sum on µ) where Λ µ (p ′ , p) is the vertex function derived in Appendix A. The external quarks take normalization factors N as well as spinor factors [7].

Chromoelectric field: µ = 4
For the interaction with the chromoelectric background field, we use the time component J 4 . To O(p 2 /m 2 ) the current in continuum QCD is where P = (p ′ + p)/2. After a short calculation with the new lattice action we find The correct (tree-level) matching is achieved if one adjusts .18) and (c E , r E ) such that m E = m 2 : At fixed m 0 a the latter prescribes a line in the (c E , r E ) plane. As before, this line becomes vertical at m 0 a = 0, fixing c E = 1 and leaving r E undetermined.
To obtain conditions on c EE , r EE , and z EE , we shall have to turn to Compton scattering in Sec. IV D.

Chromomagnetic field: µ = i
For the interaction with the chromomagnetic background field, we use the spatial components J i . To O(p 3 /m 3 ) the current in continuum QCD is After another short calculation we find where m 2 , m 3 4 , w 4 , and m 2 E have been introduced already, and (4.28) The term w X a 3 X is discussed below. Comparing Eqs.
Taken with Eqs. (4.6) and (4.7), these tuning conditions put eight constraints on the nine (nonredundant) couplings for interactions made solely out of spatial derivatives (and, hence, chromomagnetic fields). To eliminate z 6 from the right-hand side of Eq. (4.33), and to obtain conditions on r BB and z BB , we shall have to turn to Compton scattering in Sec. IV D. Equations (4.29)-(4.34) make concrete several abstract features of Sec. II. If one would like to take c 4 to be redundant in Eq. (4.7), then one cannot take r 5 to be redundant here, and similarly for z 6 and r 7 or r ′ 7 . Also, a mistuned c 5 − r 5 leads to w ′ B = 0 and a spin-dependent contribution The mismatch here is suppressed by λ 2 in the HQET counting-as expected from Table II-and by a 3 in the usual Symanzik counting.
The only undesired term in Eq. (4. 21) not yet discussed is 1 4 w X a 3 X i , where One cannot tune w X = 0. Fortunately, however, X = 0. A simple geometric proof is as follows: if, by chance, P is parallel to K, then setting P ∝ K one sees that the last two terms on the right-hand side of Eq. (4.35) vanish and the first two cancel. In the general case that P is not parallel to K, then K, P , and K × P are three linearly independent vectors. But one easily sees that K · X = P · X = (K × P ) · X = 0; (4.37) thus, X = 0. Such identities are very useful in simplifying expressions for the Compton scattering amplitude.

C. Quark-quark scattering
To match the four-quark action, Sq qqq , one must work out the quark-quark scattering amplitude. With the current J µ derived in the previous subsection, this is a relatively simple task. The main new ingredient is the improved gluon propagator. For k 2 a 2 ≪ 1, one finds [17] D µν (k) = n µ (k)D cont µν (k)n ν (k) 1 + xa 2 k 2 + O(a 4 ), (4.38) where x is the redundant coupling of the pure-gauge action, cf. Appendix C and Ref. [19]. This approximation suffices for evaluating t-channel gluon exchange. Once the bilinear action has been matched correctly, the lattice amplitude (using, say, Feynman gauge) is clearly merely where 1 and 2 label the scattered quark flavors, and both t a have uncontracted color indices. We find, therefore, that the tree-level couplings of Sq qqq are, at most, proportional to x. They can be eliminated, at the tree level, by setting x = 0, with the added benefit of simplifying the gauge action S D 2 F 2 . Note, however, that the approximation in Eq. (4.38) and, thus, Eq. (4.39), breaks down for schannel annihilation of heavy quarks. As discussed in Sec. II E, these interactions are suppressed for other reasons, so the four-quark operators needed to correct them may be neglected.

D. Compton scattering
The matching of Secs. IV A-IV C leaves four non-redundant couplings of the new action undetermined: z 6 , c EE , z EE , and z BB . To find four more matching conditions, we turn to Compton scattering. We shall proceed with the gauge-action redundant coupling x = 0.
The amplitude is whereǭ ν and ǫ µ are continuum polarization vectors, andM ab µν denotes the sum of Feynman diagrams shown in Fig. 2. The factors n ν (k ′ ) and n µ (k) appear in Eq. (4.40) to account for lattice gluons. With them one can verify that pol.
ǫ µ (k)n µ (k)n ν (k)ǭ ν (k) = −D µν (k), (4.41) as usual. We find it convenient to associate these factors with the diagrams and introduce M ab µν = n ν (k ′ )M ab µν n µ (k). Then The propagator S(q) and vertex factors Λ µ , X µν and Y µν are defined in Appendix A. The gluon propagator, to the accuracy needed, is given in Eq. (4.38), and to the same accuracy the triple-gluon vertex is (with x = 0) Note that the factors n σ (K), etc., arise naturally. Note also that K · J = k · ǫ = k ′ ·ǭ ′ = k 2 = k ′ 2 = 0, so most of the lattice artifacts in the vertex drop out. The remaining one is necessary to cancel a similar lattice artifact from the other diagrams, cf. Eqs. (B10) and (B11). We may choose the polarization vectors such thatǭ ′ 4 = ǫ 4 = 0. Then we need only focus on M mn . We have verified that M 44 is improved by (a subset of) the improvement conditions needed for A(qg → qg) calculated with these polarization vectors.
The present the results, let us introduce some notation. Write the momenta as so q = P + R and q ′ = P − R. Note that P 0 = −iP 4 = 2m 1 + · · · is larger than the other momenta, and K 0 = −iK 4 = (p ′ 2 − p 2 )/2m 2 is smaller. Next separate the diagrams according to a color decomposition, To match to continuum QCD one requires z EE = 0 (4.52) and the adjustment of (c EE , r EE ) so that m EE = m 2 . As with, say, (c E , r E ), at fixed m 0 a the latter prescribes a line in the (c EE , r EE ) plane, which becomes vertical at m 0 a = 0, fixing c EE = − 1 8 and leaving r EE undetermined. The (n, r) = (3, 0) amplitudes are where "matched" denotes terms (spelled out in Appendix B) that already match, if the conditions derived so far are applied. Equations (4.53) and (4.55) yield the new conditions

E. Matching Summary
Equations (4.6), (4.7), (4.31)-(4.34), (4.59), and (4.60) can now be combined to yield To run a numerical simulation, we would like to have as few new couplings as possible. The matching calculations verified the presence of several redundant directions. We may, therefore, take r 5 = r 7 = r ′ 7 = r BB = 0 8m 0 a(2 + m 0 a) and we also find z E = z EE = 0. , we find only six that are nonzero at tree-level matching. Moreover, once the bilinear action has been matched, and the redundant gauge coupling x = 0, the only non-zero four-quark interaction would correspond to (highly suppressed) QQ annihilation. In the next section we shall examine the size of the remaining uncertainties, to justify that this level of matching suffices.

V. ERRORS FROM TRUNCATION
In this section we give a semi-quantitative analysis of heavy-quark discretization effects with the new action. Our aim is to study the accuracy needed in matching lattice gauge theory to continuum QCD. Several elements are needed. First, we need estimates of the mismatch at short distances. This is straightforward, because the calculations of Sec. IV can be applied to work out how large the mismatch is for the unimproved action. Second, we need estimates of the longdistance effects, which is possible parametrically, by counting powers of Λ and υ. Finally, the size of discretization effects depends on the lattice spacing (obviously) so we must note the range that is tractable today and in the near future.
The error analysis is convenient using the non-relativistic description. Heavy-quark effects of operators that are related as in Eqs. (2.14) and (2.15) are lumped into one short-distance coefficient C lat i per HQET operator in Table III. In Sec. IV the short-distance coefficients are 1/2m 2 , 1/2m B , 1/4m 2 E , 1/8m 3 4 , w 4 , w B i , etc. In the corresponding continuum short-distance coefficients C cont i , these masses are replaced with a single mass m Q . To eliminate discretization effects from the kinetic energy, one should identify m Q with m 2 .
Comparison of Eqs. (2.5) and (2.6) then says that heavy-quark discretization effects take the form For example, the error from (p 2 ) 2 /8m 3 4 is See Refs. [11,12] for further details, and Ref. [31] for the application of this technique to compare several heavy-quark formalisms. We estimate the matrix elements O i using the power counting of HQET and NRQCD for heavy-light hadrons and quarkonium, respectively. The power of λ or υ is listed in Table III. The coefficient mismatches are obtained from Sec. IV, where explicit expressions show how the coefficients depend on the new couplings. In particular, when the new couplings vanish, we derive the mismatch for the Wilson and clover actions.
Explicit calculations of the mismatch at higher orders of perturbation theory are not yet available. (They would be tantamount to higher-loop matching.) Nevertheless, the asymptotic behavior remains constrained, when m Q a ≪ 1 because of the presence of the E-type operators, when m Q a ≪ 1 by heavy-quark symmetry, and when m Q a ∼ 1 because the Wilson time derivative ensures only one pole in the propagator [11]. It turns out that the most pessimistic asymptotic behavior for 1/2m B , 1/4m 2 E , etc., is the same at higher orders as in the tree level formulas in Sec. IV. It seems reasonable, therefore, to multiply the tree-level mismatch with α l s to estimate the l-loop mismatch. We use one-loop running for α s (a) starting with α s (1/11 fm) = 1/3. This yields the high end of the Brodsky-Lepage-Mackenzie coupling [32] calculated for similar quantities [33].
The resulting estimates for the mismatch of rotationally symmetric operators are shown in Fig. 3, as a function of the lattice spacing a = m 2 a/m Q , Q ∈ {c, b}. We show the relative error in mass splittings, which are of order Λ in heavy-light hadrons and of order m Q υ 2 in quarkonium. The left set of plots uses HQET power counting, for heavy-light hadrons, while the right set of plots uses NRQCD power counting, for quarkonia. The light gray or red (dark gray or blue) curves show the estimate for hadrons containing c (b) quarks. The dotted curves show the error when the corresponding correction term is omitted completely, i.e., the errors in the Wilson action. The dashed (solid) curves show the estimate of the error for tree-level (one-loop) matching. The vertical lines highlight a = 0.125 fm, 0.09 fm, 0.06 fm and 0.045 fm, corresponding to the ensembles of gauge fields with n f = 2 + 1 flavors from the MILC collaboration [34].
To drive the each contribution to heavy-quark discretization effects below 1%, we find that oneloop matching is necessary for c B , the coupling of the chromomagnetic clover term. Tree-level matching is sufficient for the chromoelectric clover coupling c E , though one-loop matching would be desirable for charmonium and charmed hadrons. The lowest plots, labeled "from 1/8m 3 4 " are for the relativistic correction terms, with couplings c 2 and z 6 . They also apply to 1/8m 3 B ′ and the related chromomagnetic couplings c 3 and z 7 . The one-loop mismatches of four-quark interactions are suppressed not only by a loop factor, but also by λ 2 or υ 2 , so they should fall below 1% too.
Similar results for operators that break rotational symmetry are shown in Fig. 4. To drive these contributions to heavy-quark discretization effects below 1%, we again find it sufficient to tune the couplings of the new action at the tree level.
There are some other noteworthy features of Figs. 3 and 4. For m Q a ≪ 1, the discretization effects vanish as a power of a, as one would deduce from the Symanzik effective field theory.   Because we identify m 2 with the mass in the C cont i , the powers of a are balanced by Λ or m Q υ, not m Q . Had we identified m 1 with the physical mass, errors of order (m Q a) n would have appeared. For m Q a ∼ 1, the tree-level curves flatten out. The error cannot grow without bound, because of the heavy-quark symmetries of the Wilson action and our improvements to it. Indeed, the curves for the b quark are usually lower than those for the c quark, which bodes well for calculations relevant to the Cabibbo-Kobayashi-Maskawa matrix. The underlying reason for the pattern is that the static approximation works better for b-flavored hadrons than for charmed hadrons. The 1/m n b contributions start out smaller, so their mismatches are also smaller. Similarly, the leading NRQCD works better for bottomonium than charmonium. The mismatches from 1/8m 3 4 and w 4 /6 deviate from the pattern, however, because NRQCD's relative suppression υ 2 bb /υ 2 cc is not as strong as HQET's (m c /m b ) 3 . Mismatches from w B i /4 and (w 4 + w ′ 4 )/4 are of order υ 4 and again follow the pattern.
In tree-level improvement, one should avoid choices where it is known that one-loop corrections from tadpole diagrams will be large [35]. Therefore, we envision following some sort of tadpole improvement. In the action, write each link matrix as u 0 [U µ /u 0 ] and absorb all but one pre-factor of u 0 into a tadpole-improved couplingc i andr i . (In several cases, it will be necessary to expand expressions such as D ilat △ ilat , △ i 2 lat , and Eq. (3.11), to eliminate any instance of U µ U † µ = 1 before inserting u 0 .) Then apply the conditions of Sec. IV toc i andr i instead of c i and r i , and take the u 0 factors in the denominator from the Monte Carlo simulation.

VI. CONCLUSIONS
In this paper we have presented the formalism and explicit calculations needed to define a new lattice action for heavy quarks. Our aim was to obtain an action whose discretization errors would be 1% at currently available lattice spacings. Combining our matching calculations, power counting, and the heavy-quark theory of discretization effects, we have argued that the proposed action should meet its target. Setting to zero the redundant couplings and those that vanish when matched at the tree level, our action can be written S = S 0 + S B + S E + S new , where The new action has six additional nonzero couplings, which depend on the couplings in S 0 + S B + S E according to Eqs. (4.73)-(4.76) and (4.79). To achieve 1% accuracy, S B must be, and S E could well be, matched at the one-loop level [36]. Another lattice action achieves similar accuracy for charmed quarks, namely the highlyimproved staggered quark (HISQ) action [37]. Our approach is computationally more demanding than HISQ. Its advantage, however, is the intriguing result that our discretization errors for bottom quarks are smaller than for charmed quarks. That means that experience with charmed hadrons and charmonium can inform analogous calculation of properties of b-flavored hadrons.
Finally, we note that there is tension between the most accurate calculation of the D s meson decay constant, f Ds [38], which uses HISQ, and experimental measurements [39]. Our action is a candidate for the charmed quark in a cross-check of the HISQ f Ds , because its discretization errors can be expected to be small enough to strengthen or dissipate the disagreement, while possessing different systematic errors.

Acknowledgments
We thank Massimo Di Pierro, Aida El-Khadra, and Paul Mackenzie for helpful conversations. Colin Morningstar provided useful correspondence on unpublished details of improved anisotropic gauge actions [23,24] The tree-level mass shell is p 4 = iE, where the energy satisfies Incoming external fermion lines receive factors u(ξ, p)N (p) or v(ξ, p)N (p), where is not easy to express in closed form. We refer the reader to two papers of Weisz for details [17] and a correction [18] for the propagator on isotropic lattices. The improved vertex is in Ref. [18]. subscript "lat" in this Appendix. One-gluon vertices need It is convenient to write out the chromomagnetic and chromoelectric cases of Eq. (A10): △ ρ ab ,µν (P, k, l) = where K = k + l. For the clover operator it is convenient to introduce C µν (k, l) = 2 cos 1 2 (k + l) µ a cos 1 2 l µ a cos 1 2 (k + l) ν a cos 1 2 k ν a − cos 1 2 k µ a cos 1 2 l ν a.
Then one has (K = k + l) The Feynman rules for one gluon are then with In the r 5 and z ′ 7 terms, Eq. (3.11) has been assumed. If instead one prefers Eq. (3.12) then replace Both choices have the same effect on Eq. (4.21). The two-gluon rules are where now p ′ = p + k + l, and s = p ′ + p = 2p + k + l;

APPENDIX B: DETAILS OF COMPTON AMPLITUDES
The parts of the Compton scattering amplitude not exhibited in Sec. IV D are shown here. First the color-symmetric contributions: 4 a 3 (w 4 + w ′ 4 ) (P n ε mri + P m ε nri )K r iΣ i − a 3 1 4 (w 4 + w ′ 4 )ε mnr K r (P m iΣ m − P n iΣ n ).
The color-antisymmetric contributions from Fig. 2(a)-(c): [P m ε nri (R − 1 2 K) r − P n ε mri (R + 1 2 K) r ]iΣ i m 2 m B − (iΣ m ε nrs − iΣ n ε mrs )R r K s + ε mnr K r iΣ · R − ε mnr R r iΣ · K 4m 2 The terms on the last line do not match, but we still must add to Eqs. (B6)-(B10) the contribution of the diagram with the three-gluon vertex [ Fig. 2(d)], which is and no M µν contribution. Here J µ is the current of Sec. IV B. The first lattice artifact cancels the last line of Eq. (B10). The second lattice artifact vanishes upon contraction with the external-gluon polarization vectors.

APPENDIX C: IMPROVED GAUGE ACTION
In this Appendix we outline how to improve the gauge action, when axis-interchange symmetry is given up. The improvement program is the same as for anisotropic lattices, which has been worked out [24] and summarized [23]. Since it has not been published, we give the main details here.
work [17][18][19] used tr[(D µ F ρν )(D µ F ρν )], but we find it easier to use tr[F µν F νρ F ρµ ]. With the Bianchi identity D µ F ρν + D ρ F νµ + D ν F µρ = 0, one can show that where ∂ denotes the omission of total derivatives that make no contribution to the action. Thus, only two of these three operators are needed. Table VI is laid out in a suggestive way: operators in the right column clearly descend from those in the left. It is a little harder to show that there are no more [24]. When parity and charge conjugation are taken into account there are 10 operators with two Ds and two Es and another 10 where the two Es are replaced with two Bs. Of these 2 × 6 may be eliminated in favor of total derivatives and others, leaving 2 × 4 = 8 of this type. Three of these may be eliminated with the Bianchi identities  Rectangles with temporal long side 1t ′ Rectangles with temporal short side 1s Spatial rectangles 2t "Parallelograms" with two temporal sides 2s Spatial "parallelograms" 3t Bent rectangles with temporal bend edge 3t ′ Bent rectangles with temporal sides, but spatial bend edge 3s Spatial bent rectangles To construct an improved gauge action, it is enough to consider the eight classes of six-link loops shown in Fig. 1, as well as plaquettes. Generalizing from Ref. [19], we label sets of unoriented loops as in Table VII. Then let where U(C) is the product of link matrices around the curve C. The gauge action is where the c i are chosen so that S D 2 F 2 ≥ 0 and so that classical continuum limit is correct. The classical continuum limit is needed not only to determine the normalization of the c i , but also to deduce which terms in the lattice action correspond to the redundant operators of the LEL. The classical continuum limit of the S i is easy to find with the procedure given in Ref. [19]. For the plaquette terms we find where a t and a s are temporal and spatial lattice spacings, respectively. Here It is convenient to express the six-link loops through S 0t and S 0s , plus further terms of order a 2 . The rectangles yield where x t , x t ′ , and x s are free parameters.