Renormalization constants of local operators for Wilson type improved fermions 2 1 0 C. Alexandrou a,b, M. Constantinou a, T. Korzec c, H. Panagopoulos a, F. Stylianou a∗ 2 n a Departmentof Physics, University of Cyprus, a PoB 20537, 1678 Nicosia, Cyprus J 4 2 b ComputationbasedScience and Technology ResearchCenter, The Cyprus Institute, 15 Kypranoros Str., 1645 Nicosia, Cyprus ] t a cInstitut fur Physik, HumboldtUniversitat zu Berlin, l - Newtonstrasse15, 12489 Berlin, Germany p e h [ (Dated: January 25, 2012) 1 Perturbative and non-perturbative results are presented on the renormalization constants of the v quarkfield and thevector,axial-vector, pseudoscalar, scalar and tensor currents. Theperturbative 5 computation,carriedoutatone-looplevelanduptosecondorderinthelatticespacing,isperformed 2 forafermionaction,whichincludestheclovertermandthetwistedmassparameteryieldingresults 0 that are applicable for unimproved Wilson fermions, as well as for improved clover and twisted 5 mass fermions. We consider ten variants of the Symanzik improved gauge action corresponding . 1 to ten different values of the plaquette coefficients. Non-perturbative results are obtained using 0 the twisted mass Wilson fermion formulation employing two degenerate dynamical quarks and the 2 tree-level Symanzik improved gluon action. The simulations are performed for pion masses in the 1 range of 480 MeV to 260 MeV and at three values of the lattice spacing, a, corresponding to : v β = 3.9, 4.05, 4.20. For each renormalization factor computed non-perturbatively we subtract its i perturbative O(a2) terms so that we eliminate part of the cut-off artifacts. The renormalization X constants are converted to MS at a scale of µ = 2 GeV. The perturbative results depend on a r large number of parameters and are made easily accessible to the reader by including them in the a distribution package of thispaper, as a Mathematica input file. PACSnumbers: 11.15.Ha,12.38.Gc,12.38.Aw,12.38.-t,14.70.Dj Keywords: LatticeQCD,Twistedmassfermions,Renormalizationconstants, Improvement ∗Electronicaddress: [email protected],[email protected],[email protected],[email protected],[email protected] 2 I. INTRODUCTION Simulations of Quantum Chromodynamics (QCD) are nowadays being carried out at almost physical parameters. Therefore studies of hadron structure within lattice QCD are beginning to yield results that can be connected to experimentmore reliably thaneverbefore. In these lattice QCDstudies one calculates matrixelements oflocaloper- atorsbetween hadronstates. Unless these operatorscorrespondto a conservedcurrentthey haveto be renormalized. Calculation of renormalization factors can be carried out using lattice perturbation theory. Although perturbation theory on the lattice is computationally more complex than in the continuum these calculations can be extended beyond one-loop order [1–3]. Various methods to improve the convergence of lattice perturbation theory have been introduced [4, 5] yielding valuable first input to the values of the renormalization constants. In this work we will use the perturbative results to improve the non-perturbative evaluation of the renormalization constants. We use the Rome-Southampton method (also known as the RI-MOM scheme) [6] to compute renormalization coefficients of arbitrary quark- antiquark operators non-perturbatively. In this approach the procedure is similar to that used in continuumperturbationtheory. Inparticular,therenormalizationconditionsaredefinedsimilarlyinperturbativeand non-perturbativecalculations. The renormalizationfactors,obtained for different values of the renormalizationscale, areevolvedperturbativelytoareferencescaleµ=2GeV.Inaddition,theyaretranslatedtoMSat2GeVusing3-loop perturbative results for the conversion factors. Since in the end one wants to make contact with phenomenological studies,whichalmostexclusivelyrefertooperatorsrenormalizedinthe MSscheme ofdimensionalregularization,one needstherenormalizationfactorsleadingfromthebareoperatorsonthelatticetotheMSoperatorsinthecontinuum. A number of lattice groups are producing results on nucleon form factors and first moments of structure functions closer to the physical regime both in terms of pion mass as well as in terms of the continuum limit [7–13]. In these lattice QCD computations one calculates hadron matrix elements of bilocal operators. In order to compare hadron matrix elements of local operators to experiment one needs to renormalize them. The aim of this paper is to calculate non-perturbatively the renormalization factors of the vector, axial-vector, scalar, pseudoscalar and tensor currents within the twisted mass formulation of Wilson lattice QCD [14]. We show that, although the lattice spacings considered in this work are smaller than 1 fm, O(a2p2) terms are non-negligible and are significantly larger than statistical errors. We therefore compute the O(a2g2)-terms perturbatively and subtract them from the non- perturbativeresults. Thissubtractionsuppresseslatticeartifactsconsiderablydependingontheoperatorunderstudy and leads to a more accurate determination of the renormalizationconstants. This approachwas applied to evaluate the renormalizationconstants for one-derivative bilinear operators in Ref. [15]. The paper is organized as follows: in Section II we give the expressions for the fermion and gluon actions we employed, and define the operators. Sections III and IV concentrate on the perturbative procedure, and the O(a2)- corrected expressions for the renormalization constants Z and Z . In Section V we provide the renormalization q O prescription of the RI′-MOM scheme, and we discuss alternative ways for its application, while in Section VI we provide all necessary formulae for the conversion to MS and the evolution to a reference scale of 2 GeV. Section VII focusesonthenon-perturbativecomputation,whereweexplainthedifferentstepsofthecalculation. Themainresults of this work are presented in Section VIII: the reader can find numerical values for the Z-factors of the fermion field and fermion operators, which are computed non-perturbatively and corrected using the perturbative O(a2) terms presented in Sections III and IV. For comparison with phenomenological and experimental results, we convert the Z-factors to the MS scheme at 2 GeV. In Section IX we give our conclusions. II. FORMULATION A. Lattice actions Our perturbative calculation makes use of the twisted mass fermion action including the usual clover (SW) term with a clover parameter that is left free. For N flavor species and using standard notation, this action reads F a3 S =− ψ¯ (x)(r−γ )U ψ (x+aµ)+ψ¯ (x+aµ)(r+γ )U ψ (x) F f µ x,x+aµ f f µ x+aµ,x f 2 x,f,µ(cid:20) (cid:21) X 4r +a4 ( +mf +iµfγ τ3)ψ¯ (x)ψ (x) a 0 0 5 f f x,f X a5 − rc ψ¯ (x)σ F (x)ψ (x), (1) SW f µν µν f 4 x,f,µ,ν X 3 wheretheWilsonparameterr ishenceforthsettor =1,f isaflavorindex,σ =[γ , γ ]/2andtheclovercoefficient µν µ ν c as well as the twisted mass parameter µf are kept as free parameters throughout. F is the standard clover SW 0 µν discretization of the gluon field tensor [16]. For the non-perturbative calculation, we consider the purely twisted mass fermion action (no clover term), which for two degenerate flavors of quarks is given by: 1 −→ −→∗ ar−→ −→∗ S =a4 χ(x)( γ (∇ + ∇ )− ∇ ∇ +m +iµ γ τ3)χ(x), (2) F 2 µ µ µ 2 µ µ 0 0 5 x X withτ3 the Paulimatrixactinginthe isospinspace,andµ the baretwistedmass. MaximallytwistedWilsonquarks 0 are obtained by setting the untwisted bare quark mass m to its critical value m , while the twisted quark mass 0 cr parameter µ is kept non-vanishing in order to give the light quarks their mass [14, 17]. The physical quantities 0 extractedfromlattice simulationsemployingmaximallytwistedWilsonquarksareautomaticallyO(a) improved[17]. In Eq. (2) the quark fields χ are in the so-called “twisted basis”. The “physical basis” is obtained for maximal twist by the simple transformation: iω iω π ψ(x)=exp γ τ3 χ(x), ψ(x)=χ(x)exp γ τ3 , ω = . (3) 5 5 2 2 2 (cid:18) (cid:19) (cid:18) (cid:19) In terms of the physical fields the action is given by: Sψ =a4 ψ(x) 1γ [−→∇ +−→∇∗]−iγ τ3 −a −→∇ −→∇∗ +m +µ ψ(x). (4) F 2 µ µ µ 5 2 µ µ cr 0 Xx (cid:18) (cid:16) (cid:17) (cid:19) One can check that this action is equivalent to the action in the twisted basis given by Eq. (2), by performing the rotations defined in Eq. (3) and identifying m =m . 0 cr For gluons we employ the Symanzik improved action, involving Wilson loops with 4 and 6 links (1×1 plaquette, 1×2 rectangle, 1×2 chair, and 1×1×1 parallelogram wrapped around an elementary 3-d cube), which is given by 2 S = c ReTr{1−U } + c ReTr{1−U } G g2 0 plaq. 1 rect. 0h pXlaq. rXect. + c ReTr{1−U } + c ReTr{1−U } . (5) 2 chair 3 paral. cXhair pXaral. i The coefficients c can in principle be chosen arbitrarily, subject to the following normalization condition, which i ensures the correct classical continuum limit of the action: c +8c +16c +8c =1. (6) 0 1 2 3 Some popular choices of values for c used in numerical simulations will be considered in this work, and are itemized i in Table I (the acronym TILW represent the Tadpole Improved Lu¨scher-Weisz action); they are normally tuned in a way as to ensure O(a2) improvement in the pure gluon sector. In our non-perturbative computation presented here we employ the tree-level Symanzik action (c = 5/3, c = −1/12, c = c = 0). Our 1-loop Feynman diagrams 0 1 2 3 do not involve pure gluon vertices, and the gluon propagator depends only on three combinations of the Symanzik parameters: C ≡c +8c +16c +8c =1, 0 0 1 2 3 C ≡c +c , (7) 1 2 3 C ≡c −c −c . 2 1 2 3 Therefore, with no loss of generality we can set c =0. 2 4 Action c c c 0 1 3 Plaquette 1.0 0 0 Symanzik 5/3 -1/12 0 TILW, βc =8.60 2.3168064 -0.151791 -0.0128098 0 TILW, βc =8.45 2.3460240 -0.154846 -0.0134070 0 TILW, βc =8.30 2.3869776 -0.159128 -0.0142442 0 TILW, βc =8.20 2.4127840 -0.161827 -0.0147710 0 TILW, βc =8.10 2.4465400 -0.165353 -0.0154645 0 TILW, βc =8.00 2.4891712 -0.169805 -0.0163414 0 Iwasaki 3.648 -0.331 0 DBW2 12.2688 -1.4086 0 TABLEI: Inputparameters c , c , c . 0 1 3 B. Definition of operators and Renormalization condition The ultra-local bi-fermion operators considered in this work are the following: ψ¯τaψ a=1,2 Oa =χ¯τaχ = (8) S (−iψ¯γ511ψ a=3 ψ¯γ τaψ a=1,2 Oa =χ¯γ τaχ = 5 (9) P 5 (−iψ¯11ψ a=3 ψ¯γ γ τ2ψ a=1 5 µ OVa =χ¯γµτaχ =−ψ¯γ5γµτ1ψ a=2 (10) ψ¯γµτ3ψ a=3 ψ¯γµτ2ψ a=1 OAa =χ¯γ5γµτaχ =−ψ¯γµτ1ψ a=2 (11) ψ¯γ5γµτ3ψ a=3 Oa =χ¯σ τaχ =ψ¯σµντaψ a=1,2 (12) T µν (−iψ¯γ5σµν11ψ a=3 ψ¯γ σ τaψ a=1,2 Oa =χ¯γ σ τaχ = 5 µν (13) Tp 5 µν (−iψ¯σµν11ψ a=3 For convenience we have included Oa even though its components are related to those of Oa. We denote the Tp T corresponding Z-factors by Za, Za, Za, Za,Za,Za . In a massless renormalization scheme the renormalization S P V A T Tp constants are defined in the chiral limit, where iso-spin symmetry is recovered. Hence Z-factors become independent of the isospin index a = 1,2,3 and we drop the a index on the Z-factors from here on. Still note that, for instance, the physicalψ¯γ τ1ψ is renormalizedwith Z while ψ¯γ τ3ψ needs Z , which differ fromeachother evenin the chiral µ A µ V limit. The renormalization constants are computed both perturbatively and non-perturbatively in the RI′-MOM scheme at different renormalization scales. We translate them to the MS scheme at µ =2 GeV using a conversion factor computed in perturbationtheory to O(g6) as describedin Section VI. The Z-factorsaredetermined by imposing the following conditions in the massless theory, i.e., at critical mass and vanishing twisted mass 1 Z = Tr (SL(p))−1S(0)(p) (14) q 12 p2=µ2 h i(cid:12) (cid:12) 1 (cid:12) Z−1Z Tr ΓL(p)Γ(0)−1(p) = 1, (15) q O 12 p2=µ2 h i(cid:12) (cid:12) (cid:12) 5 where the trace is taken over spin and color indices, µ is the renormalization scale, while SL and ΓL correspond to the perturbative or non-perturbative results, and S(0) is the tree-level result for the propagator defined as: −i γ sin(p ) S(0)(p)= ρ ρ ρ , (16) sin(p )2 Pρ ρ P while Γ(0) is the tree-level value for the fermion operators S, P, V, A, T, T′, that is Γ(0)(p)=11, γ5, γµ, γ5γµ, γ5σµν, σµν, (17) respectively. The trace is taken over spin and color indices. For alternative renormalization prescriptions the reader can refer to Ref. [15]. The choices for S(0) and Γ(0) given in Eqs. (16) - (17) are optimal, since we obtain Z = 1, Z = 1 when the q O gauge field is set to unity. Similarly, in the perturbative computation this condition leads to Z =1, and Z =1 at q O tree-level. III. CORRECTIONS TO THE FERMION PROPAGATOR The fermion propagator of the interacting theory is given by the following 2-point correlation function (Green’s function), with the various quantities computed in perturbation theory: hχa,f(x)χ¯b,g(y)i= πa d4p δab eip(x−y) S · ∞ −S−1 (p) · S n fg, (18) α β (2π)4 tree amp tree Z−πa (cid:18) nX=0(cid:0) (cid:1) (cid:19)αβ withS−1 (p)beingtheamputated,1PI2-pointfunctioninmomentumspace,computedperturbativelyuptoadesired amp order. S is the tree-level propagator for the twisted mass action and it is given by tree Stree= ip◦/+M(p)1+iµ γ5τ3, M(p)≡mf0 + a2 sin2(apµ/2), p◦/≡ γµa1sin(apµ), (19) 0 µ µ X X whereα, β areDiracindices,f , f areflavorindicesinthefundamentalrepresentationofSU(N ),anda, barecolor 1 2 F indices inthe fundamentalrepresentationofSU(N ). The dotproductruns overflavorandDirac indices. Due to the c diagonal form of the τ3 matrix, and since we are studying the case of only two degenerate quarks (up/down) we can simplify the expressionof S and omit τ3 by giving a flavor index to the twisted mass parameter, and at the same tree time we take mf →m : 0 0 1 S = , (20) tree ip/◦+M(p)+iµ(f)γ5 0 where now µ(1) = +µ for the up quark propagator and µ(2) = −µ for the down quark propagator. The 1-loop 0 0 0 0 Feynman diagrams that enter our 2-point amputated Green’s function calculation (S−1 ), are depicted in Fig. 1. amp 1 2 FIG. 1: One-loop diagrams contributing to the fermion propagator. Wavy(solid) lines represent gluons (fermions). 6 For the algebraic operations involved in evaluating Feynman diagrams, we make use of our symbolic package in Mathematica. In a nutshell, the required steps for the computation of a Feynman diagram are the following (the reader can find more details in Ref. [18]): • A preliminary expression for each diagram can be obtained by contracting the appropriate vertices, which is performed automatically once the algebraic expression of the vertices and the topology (“incidence matrix”) of the diagramarespecified. To limit the proliferationofthe algebraicexpressionswe exploitsymmetries ofthe theory,and we simplify the color dependence, Dirac matrices and tensor structures. • The O(a2) computation introduces severalcomplications, especially when isolating logarithms and Lorentz non- invariantterms, which leads to a whole family of infrareddivergent integrals. These can be reduced to a minimal set ofapproximately250’basis’integrals. Thisis achievedbyconvertingallpropagatordenominatorsto astandardform (qˆ2+M2)−1 using two kinds of subtractions, one for the fermion propagator 2 4 sin4(q /2)−4 sin2(q /2) −4m sin2(q /2) 1 1 µ µ µ µ 0 µ µ = + (21) q˜2 qˆ2 ( P (cid:16)P q˜2qˆ2 (cid:17) P ) where the denominator of the fermion propagator, q˜2, is defined as 2 1 q q˜2 = sin2(q )+ m + qˆ2 +(µ(f))2, qˆ2 =4 sin2( µ), M2 =m2+µ2, (22) µ 0 2 0 2 0 0 µ (cid:18) (cid:19) µ X X and one for the gluon propagator: D(q)=D (q)+ D(q)−D (q) =D (q)+D (q) D−1 (q)−D−1(q) D(q). (23) plaq plaq plaq plaq plaq n o n o D is the 4×4 Symanzik gluon propagator; the expression for the matrix D−1 (q)−D−1(q) , which is O(q4), is plaq independent of the gauge parameter, λ, and it can be easily obtained in clos(cid:16)ed form. Moreover,(cid:17)we have δ qˆ qˆ µν µ ν (D (q)) = −(1−λ) (24) plaq µν qˆ2 (qˆ2)2 Terms in curly brackets of Eqs. (21) and (23) are less IR divergent than their unsubtracted counterparts, by one or two powers in a. These subtractions are performed iteratively until all divergent integrals (initially depending on the fermion and the Symanzik propagator) are expressed in terms of the gluon propagator, (qˆ2 +M2)−1. The computation of the divergent integrals is performed in a non-integer number of dimensions D > 4. Ultraviolet divergences are explicitly isolated `a la Zimmermann and evaluated as in the continuum. The remainders are D- dimensional,parameter-free,zeroexternalmomentumlatticeintegralswhichcanberecastintermsofBesselfunctions, and finally expressed as sums of a pole part plus numerical constants. A small subset of the infrared divergentintegrals,shown in Appendix A, contains the most demanding ones in the list; they must be evaluated to two further orders in a, beyond the order at which an IR divergence initially sets in. As a consequence, their evaluation requires going to D > 6 dimensions. A correct way to evaluate strong divergent integrals is given in detail in a previous publication [18]. • The requirednumerical integrations of the algebraic expressions for the loop integrands are performed by highly optimized Fortranprograms;these aregeneratedby our Mathematica ‘integrator’routine. Eachintegralis expressed as a sumoverthe discrete Brillouinzone offinite lattices, with varying size L (44 ≤L4 ≤1284), and evaluatedfor all values of the Symanzik coefficients listed in Table I. • Thelastpartoftheevaluationistheextrapolationofthenumericalresultstoinfinitelatticesize. Thisprocedure entailsasystematicerror,whichis reliablyestimated,usinga complexinferencetechnique;forone-loopquantities we expect a relative error smaller than 10−7. Next, we provide the total expression for the inverse fermion propagator S−1 (p), computed up to 1-loop in pert. perturbationtheory. Hereweshouldpointoutthatfordimensionalreasons,thereisaglobalprefactor1/amultiplying our expressions for the inverse propagator, and thus, the O(a2) correction is achieved by considering all terms up to O(a3). The most general expression for the inverse propagator appears in the Mathematica file Zfactors.m (see Appendix D for notation). In the main text we provide a more compact expression, for special values of the various parameters, that is tree-level Symanzik improved gluon action, c = 0, Landau gauge (λ=0), but we keep the SW Lagrangianmass and the twisted mass parameter general. 7 ap2 a2i 6p3 S−1 =m+iµγ5+i 6p+ − (25) pert. 2 6 3M2ln[1+ p2 ] +g˜2 −13.0232725(2)i6p +m 0.5834586(2)−3ln[a2M2+a2p2]− M2 p2 ! 3M2ln[1+ p2 ] +iµγ5 8.7100834(2)−3ln[a2M2+a2p2]− M2 +a −10.69642965(5)p2 p2 ! " 6M2m2ln[1+ p2 ] 3p2 3M2 −0.8530378(3)m2−1.842911859(4)M2+ M2 + +3m2+ ln[a2M2+a2p2] p2 2 2 (cid:18) (cid:19) ! 3M2 3 3M4ln[1+ p2 ] +im6p 0.3393996(2)+ + ln[a2M2+a2p2]− M2 2p2 2 2(p2)2 ! 6M2ln[1+ p2 ] +iµmγ5 −6.68582372(4)+3ln[a2M2+a2p2]+ M2 p2 !# M4 M6 3m2p2 +a2 m 2.3547298(1)p2+2.3562747(1)m2+3.46524146(4)M2− + − 6p2 3(p2)2 2(M2+p2) " p2 11M2 p2 M6 M2ln[1+ p2 ] − +3m2+ ln[a2M2+a2p2]+ −9m2−2M2− M2 4 3 3 3(p2)2 p2 (cid:18) (cid:19) (cid:18) (cid:19) ! M4 M6 3m2p2 +iµγ5 0.70640552(8)p2+6.79538844(2)m2+1.16985307(3)M2− + − 6p2 3(p2)2 2(M2+p2) p2 2M2 p2 M2 M6 M2ln[1+ p2 ] − +3m2+ ln[a2M2+a2p2]+ −9m2− − M2 4 3 3 2 3(p2)2 p2 (cid:18) (cid:19) (cid:18) (cid:19) ! + ρ p4ρ (m+iµγ5) 1 − 2M2 + M4 − 2M6 + −1 + 2M6 M2ln[1+ Mp22] (cid:16)P (cid:17)p2 2 9p2 3(p2)2 3(p2)3 (cid:18) 3 3(p2)3(cid:19) p2 ! 9m2M2 209M4 M6 7M8 +i6p 1.1471643(7)p2−0.2145514(2)m2+1.15904388(6)M2− − − + 2p2 360p2 240(p2)2 40(p2)3 73p2 3m2 2M2 1 9m2 43M2 M4 7M6 M4ln[1+ p2 ] − + + ln[a2M2+a2p2]+ + + − − M2 360 2 3 24 2p2 72p2 12(p2)2 40(p2)3 p2 (cid:18) (cid:19) (cid:18) (cid:19) ! 67M2 M4 8M6 7M8 157 +i6p3 4.2478764(2)− + − + − ln[a2M2+a2p2] 120p2 120(p2)2 15(p2)3 30(p2)4 180 + 1 + 5M2 + 5M4 − 7M6 M4ln[1+ Mp22] + i ρ′ p4ρ′ 6p 7 + M2 + 67M4 (cid:18)2 18p2 12(p2)2 30(p2)3(cid:19) (p2)2 ! (cid:16)P p2 (cid:17) 240 48p2 72(p2)2 13M6 7M8 5 5M2 M4 7M6 M4ln[1+ p2 ] + − + − − − + M2 +O(a3,g4) 24(p2)3 12(p2)4 12 4p2 4(p2)2 12(p2)3 (p2)2 (cid:18) (cid:19) !# 8 where 6p = γρp and 6p3 = γρp3. To make the above expressions less complicated we defined m ≡ m and ρ ρ ρ ρ 0 M2 =m2+µ2. We would like to point out that the up quark propagatoris obtained by the choice µ(1) =+µ , while 0 P0 P 0 for the down quark propagatorone should choose µ(2) =−µ0. Moreover, g˜2 ≡ g126CπF2 and CF ≡ N2c2N−c1. Another byproduct of this part of the computation is the additive critical fermion mass; its general expression depends on c and the Symanzik parameters. These are terms proportional to 1/a that have been left out of SW Eq. (25) for conciseness: g˜2 1 m =− ε(1)+ε(2)c +ε(3)c2 + O(g4). (26) cr a m m SW m SW a h i The quantities ε(i) (listed in Table II) are numerical coefficients depending on the Symanzik parameters. m Action ε(m1) ε(m2) ε(m3) Plaquette -51.4347118(2) 13.73313097(5) 5.71513853(1) Symanzik -40.44324019(7) 11.94821988(5) 4.662672112(4) TILW (8.60) -34.17747288(3) 10.76516514(3) 3.998348778(3) TILW (8.45) -33.9488671(1) 10.71947605(3) 3.97345187(1) TILW (8.30) -33.6344391(1) 10.65632621(4) 3.939135834(8) TILW (8.20) -33.43979350(6) 10.61705314(7) 3.917851255(1) TILW (8.10) -33.1892274(1) 10.56629305(3) 3.890401337(1) TILW (8.00) -32.87904072(9) 10.50313393(3) 3.856345868(2) Iwasaki -26.07292275(7) 9.01533524(3) 3.1061330684(3) DBW2 -11.5127475(2) 4.9953066(1) 1.351772367(3) TABLE II: Numerical results for the coefficients ε(m1), ε(m2), ε(m3) (Eq. (26)) for different actions. The systematic errors in parentheses come from theextrapolation overfinitelattice size, L→∞. IV. CORRECTIONS TO FERMION BILINEAR OPERATORS In the context of this work we also study the perturbative O(a2) corrections to Green’s functions of local fermion operators (currents) that have the form: OΓ = χ¯fα′′′′(z)Γα′′β′′χgβ′′′′(z) . (27) Xz αX′′β′′(cid:16) (cid:17) Werestrictourselvestoforwardmatrixelements(i.e. 2-pointGreen’sfunctions,zeromomentumoperatorinsertions). The symbol Γ corresponds to the following set of products of the Euclidean Dirac matrices: 1 Γ∈{S,P,V,A,T,T′}≡{11, γ5, γµ, γ5γµ, γ5σµν, σµν}; σµν = [γµ,γν], (28) 2 forthescalarOS,pseudoscalarOP,vectorOV,axialOA,tensorOT andtensorprimeOT′ operator,respectively. The matrix elements of OT′ can be related to those of OT; this is a nontrivial check for our calculational procedure [18]. The relationship between the amputated 2-point Green functions ΛT and ΛT′ is: ΛµTν =−12 ǫµνµ′ν′ΛµT′′ν′. (29) µ′ν′ X The matrix elements of the above set of fermion bilinear operators can be obtained as: hχa,f(x)O χ¯b,g(y)i= πa d4p δabeip(x−y) S · Λpert.(p) · S fg, (30) α Γ β (2π)4 Γ Z−πa (cid:18) (cid:19)αβ 9 where Λpert.(p) is the amputated 1PI 2-pointGreen’s function of eachoperatorO , in momentum space,which upon Γ Γ contraction of indices becomes: ΛpΓert.(p)= hχfα(p) χ¯fα′(p)Γα′β′χgβ′(p) χ¯gβ(p)iamp.. (31) αX′β′ (cid:16) (cid:17) The only 1-particle irreducible Feynman diagram that enters the calculation of the above Green’s function is shown in Fig. 2. FIG. 2: One-loop diagram contributing to the bilinear operators. A wavy (solid) line represents gluons (fermions). A cross denotes the Dirac matrix Γ. In this diagram there are two fermion propagators, for which we allowed different µ-values, in order to have more general results. In other words, each of the two internal fermion lines on the left and on the right of the operator insertion(see Fig.2)canindependently representthe up ordownpropagator. For the evaluationofthe Z-factors,we keepthe twoflavorsindependent. TheamputatedGreensfunctionsoftheoperatorsO aregiveninthe Mathematica Γ fileZfactors.m. Asmentionedabove,onemaychoosethetwofermionpropagatorsofthediagramtocorrespondeither to the up or down quark, thus there are two twisted mass parameters µ(1), and µ(2). These can have any sign and the only restriction is: |µ(1)|=|µ(2)|. V. QUARK FIELD AND QUARK BILINEAR RENORMALIZATION CONSTANTS IN THE RI′-MOM SCHEME Anoperatorrenormalizationconstant(RC)canbethoughtofasthelinkbetweenitsmatrixelement,regularizedon thelattice,anditsrenormalizedcontinuumcounterpart. TheRCsoflatticeoperatorsarenecessaryingredientsinthe predictionofphysicalprobabilityamplitudesfromlatticematrixelements. Inthissectionwepresentthemultiplicative RCs, in the RI′-MOMscheme, of the quarkfield (Zpert.) andquark bilinear operators(Zpert.), obtained by using the q Γ perturbative expressions of S−1(p) and Λ (p). Γ The RI′-MOM renormalization scheme consists in imposing that the forward amputated Green function Λ (p), Γ computedinthechirallimitandatagiven(largeEuclidean)scalep2 =µ2,isequaltoitstree-levelvalue. Inpractice, the renormalizationcondition is implemented by requiring that in the chiral limit1: 1 Z−1Z V (p)| =1, V (p)≡ Tr Λ (p) · P , (33) q Γ Γ pρ=µρ Γ 4 Γ Γ h i where P are the Dirac projectors defined as follows: Γ PΓ ∈{PS,PP,PV,PA,PT,PT′}≡{11, γ5, γµ, −γ5γµ, −γ5σµν, −σµν}; (34) they arechosento obeythe relationTr[Γ·P ]≡4. The tracesarealwaystakenoverthe spinindices. The quarkfield Γ RC Z , which enters Eq. (33), is obtained by imposing, again in the chiral limit, the condition2: q 1 AsimplerversionofEq.(33)isgivenbytherelation: 1 1 Zq−1ZΓ4Tr ΛΓ(p) · ΛtΓree pρ=µρ = 4Tr ΛtΓree · ΛtΓree , (32) h i h i whereΛtΓree isthetree-levelvalueofΛΓ(p). 2 Strictlyspeaking,therenormalizationconditionofEq.(36)definesthesocalledRI′ scheme. IntheoriginalRI-MOMschemethequark fieldrenormalizationconditionreads: Zq−1−16iTr(cid:20)γµ∂Sq∂(ppµ)−1(cid:21)p2=µ2 =1. (35) ThetwoschemesdifferintheLandaugaugeattheN2LO. 10 Z−1V (p)| =1, V (p)≡−iTr a1 ργρsin(apρ) · S−1(p) . (36) q q pρ=µρ q 4 " a12P ρsin2(apρ) # We compute Z in the RI′-MOM renormalization scheme, definePd in Eq. (14), which can be Taylor expanded up q to O(a2) terms. This leads to: Z = − iTr ργρ(pρ− a62p3ρ) 1+ a2 ρp4ρ · S−1 (p) +O(a4g2,g4) q 4 p2 3 p2 1−loop "P ρ ρ Pρ ρ! # pρ=µρ P P i 6p a2 16p3 6pp4 = − Tr · S−1 (p)− − · S−1 (p) +O(a4g2,g4) . (37) 4 p2 1−loop 3 2 p2 (p2)2 1−loop " ! # pρ=µρ ThetraceistakenonlyoverspinindicesandS−1 istheinversefermionpropagatorthatwecomputedupto1-loop 1−loop and up to O(a2). We make the following definitions for convenience: p2 ≡ p2, p4 ≡ p4, 6p = γ p and ρ ρ ρ ρ ρ ρ ρ 6p3 ≡ γ p3. ρ ρ ρ P P P A very important issue is that the O(a2) terms in Z depend not only on |p|, but also on the direction of the P q renormalizationscale, p , as manifested by the presence of p4: ρ ρ ρ Vpert.(p)=−iTr 6p− a626p3 1+ a2 Pρp4ρ · S−1 (p) +O(a4g2,g4). (38) q 4 p2 3 p2 pert. " P ! # Asaconsequence,alternativerenormalizationprescriptions,involvingdifferentdirectionsofthe renormalizationscale µ =p , treat lattice artifacts differently. ρ ρ By implementing the perturbative expressions of S−1(p) and Λ (p) in Eqs. (33) and (36), we obtain the corre- Γ sponding RCs. For the following special choices (independent): tree-level Symanzik gauge action, c = 0, Landau SW gauge, and general mass m the results of the RCs under study are (for Z we have also kept µ and M = m2+µ2 q as free parameters): p Zpert. =1 (39) q 3ln[a2M2+a2p2] 3M2 3M4ln[1+ p2 ] +g˜2 −13.0232725(2)+am 0.3393996(2)+ + − M2 +a2 1.1471643(7)p2 2 2p2 2p22 ( " # h 2.1064977(2)p4 9m2M2 209M4 M6 7M8 −0.2145514(2)m2+1.15904388(6)M2+ − − − + p2 2p2 360p2 240p22 40p23 73p2 3m2 2M2 157p4 1 9m2 43M2 M4 7M6 M4ln[1+ p2 ] − + + + ln[a2M2+a2p2]+ + + − − M2 360 2 3 180p2 24 2p2 72p2 12p22 40p23 p2 (cid:18) (cid:19) (cid:18) (cid:19) p4 43M2 169M4 M6 7M8 1 35M2 M4 7M6 M4ln[1+ p2 ] + − + + − + − + + M2 p2 80p2 180p22 120p23 20p24 (cid:18)12 36p2 6p22 20p23(cid:19) p22 !#)(cid:12)(cid:12)pρ=µρ (cid:12) (cid:12) +O(a3g2,g4) (cid:12) The RCs ofthe bilinear operatorshavelengthyexpressionsandarenotshowninthe maintext; they arepresentedin Appendix B. We also include Appendix C which is related to the perturbative results appearing in our publication for RCs of one-derivative operators [15].