Nucleon strange electromagnetic form factors

The role of the strange quarks on the low-energy interactions of the proton can be probed through the strange electromagnetic form factors. Knowledge of these form factors provides essential input for parity-violating processes and contributes to the understanding of the sea quark dynamics. We determine the strange electromagnetic form factors of the nucleon within the lattice formulation of Quantum Chromodynamics using simulations that include light, strange and charm quarks in the sea all tuned to their physical mass values. We employ state-of-the-art techniques to accurately extract the form factors for values of the momentum transfer square up to 0.8~GeV$^2$. We find that both the electric and magnetic form factors are statistically non-zero. We obtain for the strange magnetic moment $\mu^s=-0.017(4)$, the strange magnetic radius $\langle r^2_M \rangle^s=-0.015(9)$~fm$^2$, and the strange charge radius $\langle r^2_E \rangle^s=-0.0048(6)$~fm$^2$.

Introduction: Strange quarks are the lightest nonvalance quarks in the nucleon and thus the most likely constituents to contribute to sea-quark dynamics. The study of strange-quark contributions to nucleon structure observables allows to uniquely identify sea-quark effects and understand virtual particle dynamics in the non-perturbative regime of Quantum Chromodynamics (QCD). A possible difference in the spatial distribution of strange and anti-strange quarks reflected by a non-zero strange electric form factor G s E (Q 2 ), and a finite strange magnetic moment µ s ≡ G s M (Q 2 =0) are key quantities describing the non-trivial composite structure of the nucleon. Parity violating electron-proton elastic scattering events probing the interference of photons and Z-bosons exchanges enable the measurement of the strange form factors and weak charge of the proton. An accurate determination of the neutral-weak vector form factor in combination with the electromagnetic form factors of the nucleon are needed in order to put constraints on new physics beyond the standard model (SM).
A number of major experiments have been measuring the parity violating form factors of the proton seeking to detect beyond the SM physics. The experimental program to study the strangeness in the proton has a long history beginning with the SAMPLE experiment [1,2] and continuing with the series of A4 experiments at the Mainz Microtron accelerator facility (MAMI) [3][4][5] and the HAPPEX [6][7][8][9] and G0 experiments [10,11] at JLab. However to date, the experimental results, although indicating non-zero values, carry large errors that make them inconclusive. This is confirmed by a recent global analysis of parity-violating elastic scattering data [12], where although a negative magnetic strange form factor is indicated, the large error still makes it consistent with zero. A review of the experimental program and results can be found in Ref. [13]. In addition, a number of phenomenological studies have been devoted to the study of the strangeness in the proton [14][15][16][17][18][19][20].
Given the current status of the experimental results, where there is no agreement even on the sign of the strange electromagnetic form factors, a first principle calculation of these key quantities is crucial. Lattice QCD provides a rigorous framework to compute nonperturbatively these quantities. However, it is only recently that efficient algorithms enable us to simulate the theory with physical values of the light quark masses and to evaluate disconnected quark loops to sufficient accuracy [21][22][23][24]. In this work, we use simulations generated with physical values of the light quark masses to evaluate accurately the strange quark loops and extract the electromagnetic strange form factors directly at the physical point.
Lattice methodology: The final results of this work are based on the analysis of an ensemble simulated with two mass degenerate light quarks, a strange and a charm quark (N f = 2+1+1) with masses tuned to their physical values [25]. We use the twisted mass formulation [26][27][28] including a clover term [29] for the simulations. The lattice volume is 64 3 × 128, m π L = 3.62, where L is the spatial lattice length and the pion mass m π = 0.1393(7) MeV and the lattice spacing a = 0.0801(4) fm determined from the nucleon mass [30]. We will refer to this ensemble as the cB211.072.64 ensemble. We use Osterwalder-Seiler strange and charm quarks with mass tuned to reproduce the Ω − baryon mass and the mass of Λ c respectively.
The nucleon matrix element of the electromagnetic operator decomposes into two CP-even form factors given by where F f 1 (q 2 ), F f 2 (q 2 ) are the Dirac and Pauli form factors with the upper index "f " indicating quark flavors. N (p, s) is the nucleon state with initial (final) momentum p(p ) and spin s(s ), with energy E N ( p) (E N ( p )) and mass m N . The momentum transfer squared is q 2 = q µ q µ where q µ = (p µ − p µ ) and u N is the nucleon spinor. Since we are interested in the strange quark contributions we take j s µ = e ss (x)γ µ s(x), where e s = −1/3. The electric and magnetic Sachs form factors can be expressed as linear combinations of the Dirac and Pauli form factors given in Euclidean space via the relations, where Q 2 = −q 2 . In order to extract the electric and magnetic strange form factors G s E and G s M in lattice QCD we need the evaluation of two-and three-point correlation functions, given by where x 0 = (0, 0) is the position at which the nucleon is created (source), x s the lattice point at which the nucleon is annihilated (sink) and x ins denotes the lattice site at which the current couples to a quark. We use projectors Γ ν taking for the unpolarized Γ 0 = 1 2 (1 + γ 0 ) and for the polarized Γ k = Γ 0 iγ 5 γ k . J N (x) = abc u a (x)[u bT (x) Cγ 5 d c (x)] is the standard interpolating field of the nucleon where u, d are the up/down quark fields and C is the charge conjugation matrix. We use Gaussian smearing [31,32] for the quark interpolating fields with APE-smeared [33] gauge links in the hopping operator in order to increase the overlap with the ground state. The parameters are optimized to yield a nucleon mean square radius of about 0.5 fm which ensures an early plateau in the two-point nucleon correlator [30].
The electromagnetic strange current j s necessarily couples to a vacuum strange quark. The contribution of the strange quark loop is given by To evaluate Eq. (7) we need to compute the sum over the diagonal spatial components of the strange quark propagator G(x; y) that would require 12 × L 3 inversions currently not feasible. Therefore, stochastic approaches [34] combined with dilution methods [35] are employed to compute such quark loops [36]. We perform a full dilution in spin and color [37] to avoid any stochastic contamination in that subspace. The elements of the propagator decay exponentially with the distance |x − y| thus dilution in space-time up to a specific distance reduces stochastic contamination entering from off-diagonal elements. This is implemented by employing the hierarchical probing technique [38] using a four dimensional coloring of distance-2 3 resulting in N Had = 512 Hadamard vectors. The coloring ensures an exact cancellation of up to 4 closest neighbors to the diagonal. We also exploit properties of the twisted mass fermions in the so-called one-end trick [39,40], which provides an increased noiseto-signal ratio [41,42]. We compute the quark loops for every time-slice t ins and construct the disconnected threepoint function at every value of t s and t ins by correlating it with 200 nucleon two-point functions using randomly distributed source positions per gauge configuration. We use 750 configurations, averaging over proton and neutron and forward and backward propagators to reach in total 600,000 measurements. We utilize the multi-grid algorithm implemented on GPUs through the QUDA software [43][44][45] to accelerate the calculation of the quark propagators.
The nucleon matrix element is extracted from an opti-mally constructed ratio [46][47][48] given by, which becomes time-independent for ∆E(t s − t ins ) 1 and ∆Et ins 1 yielding Π µ (Γ ν , p , p), where ∆E is the energy gap between the ground and the first excited state. In practice, one needs to identify the shortest time separation for which the excited states are sufficiently suppressed. We employ three methods to check the convergence to the ground state matrix element [36,49]: i) Plateau method: We use the ratio of Eq.(8) and identify a time-independent window (plateau) where we fit to a constant. We seek convergence of this value as we increase t s ; ii) Two-state fit: Takes into account the first excited state in the three-and two-point correlators entering in the ratio of Eq. (8) ; iii) Summation method: [50] Summing over the insertion time t ins in Eq. (8), excluding contact terms, leads to a linear behavior from where the slope yields the nucleon state matrix element as t s increases. For all the three methods we use the correlated χ 2 fits that take into account correlations between different insertion time slices and source-sink time separations.
There are several combinations of Π µ (Γ ν , p , p) from where G E (Q 2 ) and G M (Q 2 ) can be extracted as described in Ref. [30]. These lead to an over-constrained system of equations Π = D G, where the form factors G are extracted through a singular value decomposition of the coefficients D. For disconnected quantities we are not limited to use p = 0 since no additional inversions are needed. Given that p = 2π L n, we analyze three-point functions with | n | 2 ≤ 2, and | n| 2 ≤ 11 for G E and | n| 2 ≤ 26 for G M since the latter is in general more accurate allowing to reach higher values of the momentum.
In Fig. 1 we show the electric and magnetic form factors for two representative values of the momentum transfer squared. As can be seen, the plateau method yields results that are in agreement as t s is increased with the two-state fit and summation method. We then perform a weighted average over the plateau values in order to extract the final value. A similar behavior is observed for the other Q 2 values.
Renormalization: Since we use the local electromagnetic current we need to compute the renormalization function Z V , which is scheme and scale independent simplifying the renormalization procedure. We employ the Rome-Southampton method (RI -MOM scheme) [51] and use the momentum source approach introduced in Ref. [52] to achieve per mil statistical accuracy using O(10) configurations [53][54][55]. Discretization effects are suppressed using momenta that have the same spatial components, satisfying i p 4 i /( i p 2 i ) 2 <0.3 [56]. Furthermore, we subtract unwanted finite-a effects to O(g 2 a ∞ ) using results from lattice perturbation theory [55]. This procedure is performed using five N f =4 ensembles simulated with a range of pion masses in order to take the chiral limit. These gauge configurations are dedicatedly produced for the renormalization program using the same β value as the N f = 2 + 1 + 1 ensemble of this work. On each ensemble we compute 25 different values of the initial renormalization scale (a µ 0 ) 2 ∈ [1−7]. The dependence of Z V on the pion mass is very mild as confirmed by the fact that the coefficient of the quadratic term in m π is compatible with zero. After extrapolating to the chiral limit for each 25 values we then extrapolate to (aµ 0 ) 2 → 0 to remove any residual dependence on the RI -MOM scale. Due to the subtraction of the O(g 2 a ∞ ) artifacts, an almost constant line is obtained for Z V for (a µ 0 ) 2 ∈ [2 −7]. We obtain as our final value of Z V = 0.728(1)(4), where the first parenthesis gives the statistical error and the second the systematic coming from varying the fit window in the (aµ 0 ) 2 → 0 extrapolation. More details of the procedure can be found in Refs. [36,55,57].
Results: The results for the strange electric form factor G s E (Q 2 ) are presented in Fig. 2. The form factor is zero at Q 2 = 0 as expected and reaches a maximum at about Q 2 0.4 GeV 2 . In Fig. 3  employing the model independent z-expansion [58][59][60], using as t cut = (2m K ) 2 where m K = 486(4) MeV the kaon mass as measured in this ensemble. Since the series is expected to converge one can truncate to a k = k max and check convergence by increasing k max . For the electric form factor since G s E (0) = 0 we set a 0 = 0. We truncate the series to k max = 5 since including higher order terms has an insignificant effect on the fit. In order to stabilize the fit we use Gaussian priors for the coefficients with k > 1. Namely we set a k>1 = 0 ± w max(|a 0 |, |a 1 |), where w is the width of the Gaussian prior. We find that for w ≥ 10 the extracted values are unaffected and therefore we set w = 10 in the fit. We use a correlated χ 2 fit since the various Q 2 values are correlated. This improves the quality of the fit and the extracted values. From the fits we determine the strange magnetic moment given by the fit parameter a M 0 , which is the value of µ s ≡ G s M (0). The radii are extracted from the slope of the form factors as Q 2 → 0, namely via The extracted values are r 2 E s = −0.0048(6) fm 2 , r 2 M s = −0.015(9) fm 2 and µ s = −0.017 (4), where the error is purely statistical.
We perform the same analysis for the charm electromagnetic form factors. The electric charm form factor G c E is consistent with zero while the magnetic G c M tends to be negative albeit with large statistical errors that do not exclude zero for most Q 2 values. At the lowest available Q 2 value we find G c M (Q 2 0.051GeV 2 ) = −0.004 (2).
Comparison: Within the twisted mass formulation we have previously analyzed an N f = 2 ensemble with close to physical pion mass, namely m π = 130 MeV, lattice spacing a = 0.094(1) fm and lattice size 48 3 × 96 [24], referred to as the cA2.09.48 ensemble. However, our current analysis yields results with higher accuracy. Besides these two analyses, currently there are no other lattice QCD calculations of these form factors directly at the physical pion mass. The fact that we achieved the current accuracy is due to our improved methods for computing the quark loops leading to about four times more accurate results. Three other groups have computed the strange form factors with the χQCD collaboration including an ensemble with close to physical pion mass. The analysis was performed using a mixed setup with N f = 2 + 1 gauge configurations produced using domain wall fermions and overlap fermions used for the evaluation of nucleon two-and three-point correlators. The four ensembles spanned pion masses m π ∈ [139 − 330] MeV. Their final values are extracted using a chiral extrapolation since their results at the physical point alone are not accurate [22]. The other two groups used simulations with heavier than physical pions: The LHPC collaboration analyzed one ensemble of N f = 2+1 clover-improved Wilson fermions with m π = 317 MeV [21] and used an interpolation to estimate the value at the physical point; The third group [61] analyzed several CLS ensembles of N f = 2 + 1 O(a)-improved Wilson fermions with pion masses m π ∈ [200 − 360] MeV and performed a chiral extrapolation to extract the value at the physical point.
In Fig. 4, we show a comparison of the magnetic moment and radii using the two twisted mass ensembles with the corresponding results from the aforementioned groups. As can be seen, there is an overall agreement. Our very precise values for the electric radius and magnetic moment clearly confirm a non-zero value for both. The agreement among lattice QCD results using ensembles of different values of the lattice spacings and volumes also indicates that cut-off and finite volume effects are small. This allows us to make a comparison of our results obtained using the N f = 2 and the N f = 2 + 1 + 1 ensembles to check for unquenching effects of the strange quark. The current statistical accuracy reveals no such effects. Results using the N f = 2 cA2.09.48 twisted mass ensemble are shown by the blue filled square [24]. We denote with open symbols results that include ensembles with larger than physical pion masses to extract the value at the physical point. Results from the χQCD [22] collaboration (purple upper triangles), Ref. [61] (black right triangles) and from the LHPC [21](green circles). The inner error bars indicate the statistical while the outer the total which includes systematic errors.
Conclusions: A high precision calculation of the strange nucleon electromagnetic form factors is obtained using ensembles simulated with physical pion mass. Using the model independent z-expansion to fit the form factors we obtain the following values for the radii and magnetic moment clearly excluding a zero value for all three quantities. This is a significant finding given the status of experimen- 1 GeV extracted in this work. The ellipses indicate 95% confidence level. The green ellipse is from Ref. [62], the orange from Ref. [12], the black from Ref. [63] and the blue from Ref. [64]. tal searches where the results are inconclusive. For example the SAMPLE experimental data [1] finds a strange magnetic moment of µ s = 0.37 ± 0.20 ± 0.26 ± 0.15 that is positive but also compatible with zero. More recently, the HAPPEX collaboration finds at Q 2 ∼ 0.62 GeV 2 [9] a negative value for G s M = −0.070 ± 0.067, which again does not exclude zero. The G0 collaboration reported an upper bound of 10% on the strange quark contributions as compared to the total nucleon electromagnetic form factors [11]. The A4 experiment [5], also reported results consistent with zero strangeness, namely G s E = 0.050 ± 0.038 ± 0.019 and G s M = 0.14 ± 0.11 ± 0.11 at Q 2 = 0.22 GeV 2 . The Q-weak experiment [65][66][67] is aiming to measure the weak charge of the proton to unprecedented accuracy to set limits on new physics. The strange electromagnetic form factors are a crucial input that can aid the interpretation of the experimental results. In addition, the MESA [68] facility at Mainz targets very low Q 2 in order to improve the determination of the strange electric form factor, making these results of high relevance. In Fig. 5 we show the impact of the determination of the strange form factors on the experimental measurements at a given value of Q 2 = 0.1 GeV 2 . Our values provide a stringent constrain on experimental searches.
Horizon 2020 research and innovation program of the European Commission under the Marie Sk lodowska-Curie grant agreement No 642069. S.B. is supported by the Marie Sk lodowska-Curie grant agreement No 642069 of the European commission and from the project COM-PLEMENTARY/0916/0015 funded by the Cyprus Research Promotion Foundation. The project used resources of the SuperMUC at Leibniz Supercomputing Centre under the Gauss Centre for Supercomputing e.V. project pr74yo and of Piz Daint at the Centro Svizzero di Calcolo Scientifico under the project s702. It also used XSEDE computational resources supported by National Science Foundation grant number TG-PHY170022 as well as on the Jureca system of the research center in Jülich, under the project ECY00.