Strange and Charm Quark Spins from Anomalous Ward Identity
We present a calculation of the strange and charm quark contributions to the nucleon spin from the anomalous Ward identity (AWI). It is performed with overlap valence quarks on 2+1-flavor domain-wall fermion gauge configurations on a lattice with the light sea mass at MeV. To satisfy the AWI, the overlap fermion for the pseudoscalar density and the overlap Dirac operator for the topological density, which do not have multiplicative renormalization, are used to normalize the form factor of the local axial-vector current at finite . For the charm quark, we find that the negative pseudoscalar term almost cancels the positive topological term. For the strange quark, the pseudoscalar term is less negative than that of the charm. By imposing the AWI, the strange at is obtained by a global fit of the pseudoscalar and the topological form factors, together with and the induced pseudoscalar form factor at finite . The chiral extrapolation to the physical pion mass gives .
pacs:12.38.Gc, 14.20.Dh, 11.30.Hv, 14.65.Dw
The quark spin content of the nucleon was found to be much smaller than that expected from the quark model by the polarized deep inelastic lepton-nucleon scattering experiments and the recent global analysis reveals that the total quark spin contributes only to the proton spin deFlorian:2009vb ().
In an attempt to understand the smallness of the quark spin contribution from first principles, several lattice QCD calculations Dong:1995rx (); Fukugita:1994fh () have been carried out since 1995 with the quenched approximation or with heavy dynamical fermions Gusken:1999as (). The most challenging part of the lattice calculation is that of the disconnected insertion of the nucleon three-point functions due to the quark loops. Recently, the strange quark spin has been calculated with the axial-vector current on light dynamical fermion configurations QCDSF:2011aa (); Babich:2010at (); Engelhardt:2012gd (); Abdel-Rehim:2013wlz (); Chambers:2015bka () and it is found to be in the range from to . This is about 4 to 5 times smaller in magnitude than that from a global fit of DIS which gives deFlorian:2009vb () and a most recent analysis Leader:2014uua () including the JLab CLAS high precision data which finds it to be Stamenov+Leader2015 ().
Such a discrepancy between the global fit of experiments and the lattice calculation of the quark spin from the axial-vector current is unsettling. It was emphasized some time ago that it is essential that a lattice calculation of the flavor-singlet axial-vector current be able to accommodate the triangle anomaly Karsten:1980wd (); Lagae:1994bv (). It was specifically suggested Karsten:1980wd () to calculate the triangle anomaly from the VVA vertex and take it as the normalization condition for the axial-vector current in order to determine the normalization factor on the lattice. To address the discrepancy of the strange quark spin, we shall use the anomalous Ward identity (AWI) to provide the normalization and renormalization conditions to calculate the strange and charm quark spins in this work.
The anomalous Ward identity (AWI) is usually referred to the flavor-singlet axial current in the flavor basis where there is a anomaly term in the divergence of . However, the flavor is a global symmetry, AWI is satisfied for each flavor in the flavor basis through linear combinations of the flavor-octet axial current and isovector axial current . For the case of the strange quark, its AWI can be obtained from the AWI for the and the WI for (N.B. there is no anomaly term in the WI for ) through the combination . Alternatively, the AWI can be derived for the strange by considering the infinitesimal local chiral transformation where with the matrix in flavor space , where is the 8th generator, gives a chiral transformation only for the strange.
For the overlap fermion Neuberger:1997fp () which has chiral symmetry on the lattice via the Ginsparg-Wilson relation, the conserved flavor-singlet axial current is derived Hasenfratz:2002rp (). Following the derivation with the above definition for the matrix for the chiral transformation, it is straightforward to show that the following identity is satisfied for the strange axial-vector current.
where is the forward derivative. The expression of the conserved current for the strange quark is given in Ref. Hasenfratz:2002rp () which involves a non-local kernel which is more involved to implement numerically than the local current. In this work we shall replace it with the local current where is the massless overlap operator which is exponentially local with a fall-off of one lattice spacing Draper:2005mh (). The topological charge is derived in the Jacobian factor from the fermion determinant under the chiral transformation which is equal to in the continuum Kikukawa:1998pd (), i.e.
where is the trace over both spin and color, while is the trace over color. For the strange quark spin, we shall consider the in Eq. (1) to be the nucleon propagator, i.e.
where is the commonly used proton interpolation operator which involves two and one fields
where the Latin letters denotes the color index and the Greek letters denotes the Dirac index and for the Pauli-Sakurai representation that we adopt for the matrices. In this case, the first term in Eq. (1) vanishes, since does not involve strange quarks and hence no dependence.
Following the standard calculation of off-forward nucleon matrix element Liu:1994dr (); Deka:2013zha (), one considers the appropriate combination of the three-point function with the momentum projection of the current and the two-point functions to remove the kinematic dependence and take the time separation between the nucleon source and the current insertion, likewise between the nucleon sink and the current insertion, one arrives at the following unrenormalized AWI in nucleon matrix element for the strange quark
where is the nucleon state with momentum and spin . As we mentioned above that we shall replace the conserved axial-vector current with the local one . To compensate for the replacement, a normalization factor is introduced to make the AWI satisfied at finite cutoff. This is the only normalization factor needed since the pseudoscalar density and the topological charge are the same as those in Eq. (1) (N.B. In the case of the disconnected insertion for the strange quark, the pseudoscalar density contributes through the quark loop. In this case, the takes the form ). This lattice normalization factor is analogous to that introduced to make the chiral Ward identity satisfied for the local non-singlet axial-vector current. In the literature, it is usually denoted as which is actually a finite renormalization with no logarithmic scale dependence. Following Ref. Luscher:1996jn (), we shall call it lattice normalization. Unlike the vector current and non-singlet axial current, the flavor-singlet axial current has, in addition, a renormalization with anomalous dimension. We thus consider the renormalization on top of normalization as is done for the energy-momentum tensor in Ref. Deka:2013zha (). We will discuss the renormalization after we define the strange quark spin first.
The normalized strange quark spin in the nucleon , where is the bare forward matrix element from the local axial-vector current
can be obtained by evaluating the right-hand-side of the AWI in Eq. (1) between the nucleon states in the forward limit, i.e.
where and are form factors at as defined in Eq. (7). The normalized charm spin is similarly defined. In this case, one can, in principle, calculate and at finite and extrapolate them to the limit and this approach has been studied before Mandula:1990ce (); Altmeyer:1992nt (); however, the pseudoscalar density term was not included. Despite the fact that there is no massless pseudoscalar pole in the flavor-singlet case, it is shown that the contribution of the pseudoscalar density does not vanish at the massless limit Liu:1991ni (); Liu:1995kb (). Furthermore, there is a pion pole in the disconnected insertion of to cancel that in the connected insertion to lead to and poles Liu:1991ni (); Liu:1995kb (). Thus, the and form factors at small of the order of are essential for a reliable extrapolation. Since the smallest is larger than on the lattice we work on, a naive extrapolation of in Eq. (7) may lead to a wrong result. To alleviate this concern, we shall consider instead, in this work, matching the and the induced pseudoscalar form factor from the left side of Eq. (Strange and Charm Quark Spins from Anomalous Ward Identity) and and from the right side at finite to determine the normalization constant as will be discussed later.
As far as renormalization is concerned, we note that in the continuum calculation Espriu:1982bw (), the renormalization constants of the quark mass and the pseudoscalar density cancel, i.e. , and the renormalized topological charge term has a mixing with the divergence of the axial current at one-loop in the form where is the flavor-singlet axial current and with one of the coming from the definition of the topological charge. On the other hand, the renormalization of the divergence of axial-vector current occurs at the two-loop level involving a quark loop in the disconnected insertion which gives Espriu:1982bw ()the divergence of the renormalized strange axial-vector
In the present work, we adopt the overlap fermion for the lattice calculation where and there is no multiplicative renormalization of the topological charge defined by the overlap operator in Eq. (2). After two-loop matching from the lattice to the scheme, the renormalized and normalized AWI equation at the scale is therefore
where is the anomalous dimension. We see that, modulo the possible different finite terms and in the renormalization of and the topological charge Capitani:2002mp (), the anomalous dimension term on the l.h.s is the same as that on the r.h.s. Espriu:1982bw (). Thus, the two-loop renormalized AWI is the same as the unrenormalized AWI in Eq. (Strange and Charm Quark Spins from Anomalous Ward Identity).
Two loop renormalization on the lattice is quite involved, we plan to carry out the calculation of the lattice matching to the scheme non-perturbatively as is recently done in Ref. Chambers:2015bka (). For the present work, we shall give an estimate of the renormalization correction. From the left side of Eq. (Strange and Charm Quark Spins from Anomalous Ward Identity), one finds the renormalized
where is the normalized and
where is the flavor-singlet .
To estimate the size of for renormalization and matching to the scheme at GeV, we note that for the Iwasaki gauge action for DWF configurations, the lattice spacing GeV, and we assume to be . In this case, . Taking the experimental value of from experiments deFlorian:2009vb (), we obtain . We shall take this as a part of the systematic error.
We use the overlap fermion for the valence quarks in the nucleon propagator as well as for the quark loops on flavor domain-wall fermion (DWF) configurations on a lattice with the light sea quark mass corresponding to a pion mass at 330 MeV Aoki:2010dy (). Both DWF and overlap fermions have good chiral symmetry and it is shown that , which is a measure of mismatch in mixed action, is very small Lujan:2012wg () and its effects on the nucleon properties have not been found to be discernible Gong:2013vja (). Since the discretization errors are found to be small in the study of the charmonium spectrum and Yang:2014sea (), this allows us to compute the spin for the charm quark on this lattice. In addition to the advantage in normalization as mentioned above, the zero mode contributions to in the disconnected insertion (DI) and in Eq. (7), which are finite volume artifacts, cancel when the overlap operators are used for both of them.
We adopt the sum method Maiani:1987by (); Deka:2008xr () where the ratio is taken and the insertion time of the quark loop and the topological charge is summed between and where is the nucleon source/sink time. As a result, the ratio R(), where , is linearly dependent on and the slope is the matrix element of the spin content from or ,
from which we can obtain and as functions of the momentum transfer squared .
As explained in detail Gong:2013vja (); Li:2010pw (), we adopt the -noise grid smeared source, with support on some uniformly spaced smeared grid points on a time slice, and low-mode substitution (LMS) which improves the signal-to-noise ratio substantially, given the same computer resources. For the lattice, we place two smeared sources in each spatial direction, each with a Gaussian smearing radius of lattice spacing, and have seen a gain of times of statistics in the effective nucleon mass as compared to that of one smeared source. In view of the fact that the useful time window for the nucleon correlator is less than 14 and we have slices in time, we put two grid sources at and 32 simultaneously to gain more statistics from one inversion. Thus, our grid has the pattern of with two smeared grid sources in each of the space and time directions.
Since both the strange and charm are from the disconnected insertion (DI), the calculation involves the product of the nucleon propagator and the quark loop. For the quark loop, we employ the low mode average (LMA) algorithm which entails an exact loop calculation for the low eigenmodes of the massive overlap fermion over all space time points on the lattice. On the other hand, the high modes of the quark loops are estimated with 4-D noise grid sources on grids and diluted for time slices and even-odd sites for a total of 4 inversions, each with one noise.
The AWI splits the divergence of the axial current into two parts, i.e. and , and the two parts reveal different aspects of the physics contribution. The pseudoscalar part is low-mode dominated for light quarks, where the first 200 pairs of overlap eigenvectors contribute more than 90% of the vacuum value for the very light quarks and for the strange Gong:2013vja (). The overlap Dirac operator in the definition of the topological term in Eq. (2) is exponentially local with an exponential falloff of one lattice spacing Draper:2005mh (). Thus, the anomaly part, being local, captures the high-mode contribution of the divergence of the axial-vector current.
We first show the ratio in Eq. (12) for the charm quark as a function of for the case with lowest momentum transfer, i.e. GeV (corresponding to ) in the upper panel of Fig. 1. The contributions from the low modes and high modes for at this , which are coded in the slopes, are shown separately. They are from the case where the valence quark in the nucleon and that of the light sea have the same mass which correspond to MeV. It is clear from the upper panel of Fig. 1 that low modes dominate the contributions. Even though the low modes contribute only in the charm quark loop itself Gong:2013vja (), they become dominant when correlated with the nucleon. On the other hand, the from the slope at this is large and positive. The errors for and are 6% and 4% respectively.
In the lower panel of Fig. 1, we give the results for the charm quark () which is determined from a global analysis of the charm mass Yang:2014sea (). The pseudoscalar density term (red points) and the topological charge density term (black squares) are plotted as a function of . We see that the pseudoscalar contribution is large, due to the large charm mass, and negative while the anomaly is large and positive. The lines are fits with a dipole form just to guide the eye. When they are added together (blue triangles in the figure), they are very close to zero, with small statistical errors, over the whole range of . Thus, when extrapolated to with a constant, we obtain at the unitary point. When extrapolated to the physical pion mass, , which shows that the charm hardly contributes anything, if at all, to the proton spin due to the cancellation between the pseudoscalar term and the topological term. It is known Franz:2000ee () that the leading term in the heavy quark expansion of the quark loop of the pseudoscalar density, i.e. , is the topological charge which cancels the contribution from the topological term in the AWI. To the extent that the charm is heavy enough such that the correction is small, the present results of cancellation can be taken as a cross check of the validity of our numerical estimate of the DI calculation of the quark loop as well as the anomaly contribution. The mixing for the heavy quark loops from the other favors are also highly suppressed and negligible at the present stage.
Next, we consider the case with the strange quark () for this lattice, which is again determined from the global fit for the strange quark mass based on fitting of and Yang:2014sea (). Similarly to Fig. 1, and are plotted in Fig. 2 for the unitary case where the valence quarks in the nucleon and the light sea quarks have the same mass at MeV. We see in the upper panel that the low modes completely dominate the contribution as in the case of charm. The anomaly is the same for all flavors. In the lower panel, it is shown that the contribution from is only slightly smaller than that of the charm. This is due to the fact that even though the strange quark mass is about 12.5 times smaller than that of the charm Yang:2014sea (), its pseudoscalar matrix element is much larger than that of the charm. Since the anomaly is the same for the strange and the charm, the sum of and , shown in the lower panel, is slightly positive in the range of as plotted.
Since our smallest is larger than which should be present as the pion pole on the right hand side of the DI of AWI form factors to cancel that in the CI Liu:1991ni (); Liu:1995kb (), taking the limit in Eq. (7) can lead to large systematic error. In view of this, we calculated the unnormalized and the induced pseudoscalar form factor with the 3-point to 2-point correlator ratio Deka:2013zha ()
where and denote the directions of the axial current and the nucleon polarization. Here and are normalized form factors. Sandwiching the AWI between the nucleon states with finite momentum transfer, one obtains
With 18 data points for for different and 6 data points for and for 6 different , we fit Eqs. (13) and (14) to obtain (including ), , and . Since it is a global fit with all the data included, this method does not require modeling the behavior with any assumed functional form.
The results for normalized , are plotted in Fig. 3 as a function of . Also plotted is which is compared to from the AWI in Eq. (14). We see that the agreement is good for the range of except for the last point at where there is a two-sigma difference.
From the fit, we obtain
and at the unitary point where MeV. and have been calculated this way for several valence quark masses in the nucleon while keeping the quark loop at the strange quark point. The valence mass dependence of is plotted in Fig. 4. We see that is larger than 1, and becomes larger as the valence quark mass decreases.
The chiral behavior of is plotted in Fig. 5 as a function of according to the valence quark mass. We see that the results are fairly linear in . Thus we fit it linearly in with the form where is the physical pion mass and obtain at the physical pion mass. This is shown in Fig. 5. The uncertainty estimated through the variance from several different fits by adding a term, a term, or a term to the chiral extrapolation formula gives a systematic error of 0.0013.
In this work, we adopted the sum method to extract the matrix elements. To assess the excited state contamination, we shall use the combined two-state fit with the sum method used in the calculation of the and strange sigma terms Yang:2015uis (), strange magnetic moment Sufian:2016pex (), and glue spin Yang:2016plb () for comparison for a few cases. We first plot in Fig. 6 the un-summed ratios in Eq. (12) for at the smallest as a function of for time separations between the source and the sink. A combined two-state and sum method fit with these data a value of 0.035(3) which is consistent within one sigma with the slope from the sum method which is 0.033(4).
Similarly, we have done the comparison for . Plotted in Fig. 7 is the summed ratio of 3-pt-to-2-pt correlators as a unction of for the calculation of which is extracted from the slope as is from Eq. (12). At unitary point, we obtain . Also plotted in Fig. 8 are the un-summed ratios for as a function of for time separations between the source and the sink. A combined two-state and sum method fit with these data a value of -0.030(5). While their errors touch, this is larger than that from the sum method fit. We shall take this 10% difference as a systematic error of the present work.
The total systematic error contains the renormalization uncertainty , the uncertainty of the chiral extrapolation of 0.0013, and uncertainty due to the excited state contamination of the sum method of 0.0040. We sum them up quadratically and obtain an overall systematic error of 0.0078.
We list our result in Fig. 9 together with other recent lattice results in comparison with the global fit of the DIS data Leader:2014uua (); Stamenov+Leader2015 (). The blue triangles are lattice calculations of the axial vector current matrix element and the red circle is from the present work based on the anomalous Ward identity.
We see that our result, although still more than two sigmas smaller than the recent analysis of the DIS data which finds the strange spin to be -0.106(23) Stamenov+Leader2015 (), is somewhat larger in magnitude than the other direct calculations of the axial-vector current QCDSF:2011aa (); Babich:2010at (); Engelhardt:2012gd (); Abdel-Rehim:2013wlz (); Chambers:2015bka (). This is mainly due to the fact that the normalization factor , which is required to have the AWI satisfied in our calculation, is larger than that for the isovector axial-vector current which is in our case. Presumably, a similarly larger exists for the other calculations using axial-vector currents which do not satisfy the AWI, but has not been taken into account.
In summary, we have carried out a calculation of the strange and charm quark spin contributions to the spin of the nucleon with the help of the anomalous axial Ward identity. This is done with the overlap fermion for the nucleon and the quark loop on flavor DWF configurations on a lattice with light sea quarks corresponding to MeV. Since the overlap fermion is used for the pseudoscalar term and the overlap Dirac operator is used for the local topological term, the normalized AWI also holds for the renormalized AWI to two loop order. For the charm quark, we find that the and the anomaly contributions almost cancel. For the strange quark, the term is somewhat smaller than that of the charm. Fitting the AWI at finite and the and form factors, we obtain the normalized . The normalization factor for the local axial-vector current is found to be larger than that for the isovector axial-vector current, which implies that it is affected by a large cutoff effect presumably due to the triangle anomaly. This will be clarified by future work using the conserved axial-vector current Hasenfratz:1998ri () for the overlap fermion. After chiral extrapolation to the physical pion mass for the nucleon, we obtain which is consistent with zero, and which is smaller in magnitude than that from the latest analysis of DIS data Leader:2014uua (); Stamenov+Leader2015 () by more than two sigmas. We will check to see if this can be understood with lattices at the physical point. In this work, we have identified the source for the negative spin contribution in the disconnected insertion of the light quarks as due to the large and negative contribution which overcomes the positive anomaly contribution to give an overall negative . This is likely the cause for the smallness of the net quark spin in the nucleon. We will confirm this later with results on the and quarks from both the disconnected and connected insertions.
We thank RBC and UKQCD for sharing the DWF gauge configurations that we used in the present work. This work is supported in part by the National Science Foundation of China (NSFC) under the project No. 11405178, the Youth Innovation Promotion Association of CAS (2015013), and the U.S. DOE Grant DE-SC0013065. A.A. is supported in part by the National Science Foundation CAREER grant PHY-1151648. This research used resources of the Oak Ridge Leadership Computing Facility at the Oak Ridge National Laboratory, which is supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC05-00OR22725. This work also used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation grant number ACI-1053575.
- (1) D. de Florian, R. Sassot, M. Stratmann and W. Vogelsang, Phys. Rev. D 80, 034030 (2009) [arXiv:0904.3821 [hep-ph]].
- (2) S. J. Dong, J. -F. Lagae, K. F. Liu, Phys. Rev. Lett. 75, 2096-2099 (1995) [hep-ph/9502334].
- (3) M. Fukugita, Y. Kuramashi, M. Okawa and A. Ukawa, Phys. Rev. Lett. 75, 2092 (1995) [hep-lat/9501010].
- (4) S. Gusken et al. [TXL Collaboration], Phys. Rev. D 59, 114502 (1999).
- (5) G. S. Bali et al. [QCDSF Collaboration], Phys. Rev. Lett. 108, 222001 (2012) [arXiv:1112.3354 [hep-lat]].
- (6) R. Babich, R. C. Brower, M. A. Clark, G. T. Fleming, J. C. Osborn, C. Rebbi and D. Schaich, Phys. Rev. D 85, 054510 (2012) [arXiv:1012.0562 [hep-lat]].
- (7) M. Engelhardt, Phys. Rev. D 86, 114510 (2012) [arXiv:1210.0025 [hep-lat]].
- (8) A. Abdel-Rehim, C. Alexandrou, M. Constantinou, V. Drach, K. Hadjiyiannakou, K. Jansen, G. Koutsou and A. Vaquero, Phys. Rev. D 89, no. 3, 034501 (2014) [arXiv:1310.6339 [hep-lat]].
- (9) A. J. Chambers et al., [arXiv:1508.06856 [hep-lat]].
- (10) E. Leader, A. V. Sidorov and D. B. Stamenov, Phys. Rev. D 91, no. 5, 054017 (2015) [arXiv:1410.1657 [hep-ph]].
- (11) D. Stamenov and E. Leader, private communication.
- (12) L. H. Karsten and J. Smit, Nucl. Phys. B 183, 103 (1981).
- (13) J. F. Lagae and K. F. Liu, Phys. Rev. D 52, 4042 (1995) [hep-lat/9501007].
- (14) H. Neuberger, Phys. Lett. B 417, 141 (1998) [hep-lat/9707022].
- (15) P. Hasenfratz, S. Hauswirth, T. Jorg, F. Niedermayer and K. Holland, Nucl. Phys. B 643, 280 (2002), doi:10.1016/S0550-3213(02)00742-3, [hep-lat/0205010].
- (16) T. Draper et al., PoS LAT 2005, 120 (2006), [hep-lat/0510075].
- (17) Y. Kikukawa and A. Yamada, Phys. Lett. B 448, 265 (1999) [hep-lat/9806013]; D. H. Adams, Annals Phys. 296, 131 (2002) [hep-lat/9812003]; K. Fujikawa, Nucl. Phys. B 546, 480 (1999) [hep-th/9811235]; H. Suzuki, Prog. Theor. Phys. 102, 141 (1999) [hep-th/9812019].
- (18) K. F. Liu, S. J. Dong, T. Draper and W. Wilcox, Phys. Rev. Lett. 74, 2172 (1995), doi:10.1103/PhysRevLett.74.2172, [hep-lat/9406007].
- (19) M. Luscher, S. Sint, R. Sommer and H. Wittig, Nucl. Phys. B 491, 344 (1997) doi:10.1016/S0550-3213(97)00087-4 [hep-lat/9611015].
- (20) The typical form for the ratio in Eq. (13) can be found in M. Deka, T. Doi, Y.B. Yang, B. Chakraborty, S.J. Dong, T. Draper, M. Glatzmaier, M. Gong, H. W. Lin, K. F. Liu, D. Mankame, N. Mathur, and T. Streuer, Phys. Rev. D 91, no. 1, 014505 (2015) [arXiv:1312.4816 [hep-lat]].
- (21) D. Espriu and R. Tarrach, Z. Phys. C 16, 77 (1982).
- (22) For lattice perturbation theory, see for example, S. Capitani, Phys. Rept. 382, 113 (2003) [hep-lat/0211036].
- (23) J. E. Mandula, Phys. Rev. Lett. 65, 1403 (1990).
- (24) R. Altmeyer, M. Gockeler, R. Horsley, E. Laermann and G. Schierholz, Phys. Rev. D 49, R3087 (1994).
- (25) K. F. Liu, Phys. Lett. B 281, 141 (1992).
- (26) K. F. Liu, hep-lat/9510046.
- (27) Y. Aoki et al. [RBC and UKQCD Collaborations], Phys. Rev. D 83, 074508 (2011) [arXiv:1011.0892 [hep-lat]].
- (28) M. Lujan et al., Phys. Rev. D 86, 014501 (2012) [arXiv:1204.6256 [hep-lat]].
- (29) M. Gong et al. [QCD Collaboration], Phys. Rev. D 88, no. 1, 014503 (2013) [arXiv:1304.1194 [hep-ph]].
- (30) Y. B. Yang, Y. Chen, A. Alexandru, S. J. Dong, T. Draper, M. Gong, F. X. Lee and A. Li, K.F. Liu, Z. Liu, and M. Lujan Phys. Rev. D 92, no. 3, 034517 (2015) [arXiv:1410.3343 [hep-lat]].
- (31) L. Maiani, G. Martinelli, M. L. Paciello and B. Taglienti, Nucl. Phys. B 293, 420 (1987).
- (32) M. Deka, T. Streuer, T. Doi, S. J. Dong, T. Draper, K. F. Liu, N. Mathur and A. W. Thomas, Phys. Rev. D 79, 094502 (2009) [arXiv:0811.1779 [hep-ph]].
- (33) A. Li et al. [QCD Collaboration], Phys. Rev. D 82, 114501 (2010) [arXiv:1005.5424 [hep-lat]].
- (34) M. Franz, M. V. Polyakov and K. Goeke, Phys. Rev. D 62, 074024 (2000) [hep-ph/0002240].
- (35) Y. B. Yang et al. [xQCD Collaboration], Phys. Rev. D 94, no. 5, 054503 (2016), doi:10.1103/PhysRevD.94.054503, [arXiv:1511.09089 [hep-lat]].
- (36) R. S. Sufian, Y. B. Yang, A. Alexandru, T. Draper, K. F. Liu and J. Liang, Phys. Rev. Lett. 118, no. 4, 042001 (2017), doi:10.1103/PhysRevLett.118.042001, [arXiv:1606.07075 [hep-ph]].
- (37) Y. B. Yang, R. S. Sufian, A. Alexandru, T. Draper, M. J. Glatzmaier, K. F. Liu and Y. Zhao, arXiv:1609.05937 [hep-ph].
- (38) P. Hasenfratz, V. Laliena and F. Niedermayer, Phys. Lett. B 427, 125 (1998) doi:10.1016/S0370-2693(98)00315-3 [hep-lat/9801021].