The pion vector form factor from Lattice QCD at the physical point

We present an investigation of the electromagnetic pion form factor, $F_\pi(Q^2)$, at small values of the four-momentum transfer $Q^2$ ($\lesssim 0.25$ GeV$^2$), based on the gauge configurations generated by European Twisted Mass Collaboration with $N_f = 2$ twisted-mass quarks at maximal twist including a clover term. Momentum is injected using non-periodic boundary conditions and the calculations are carried out at a fixed lattice spacing ($a \simeq 0.09$ fm) and with pion masses equal to its physical value, 240 MeV and 340 MeV. Our data are successfully analyzed using Chiral Perturbation Theory at next-to-leading order in the light-quark mass. For each pion mass two different lattice volumes are used to take care of finite size effects. Our final result for the squared charge radius is $\langle r^2 \rangle_\pi = 0.443~(29)$ fm$^2$, where the error includes several sources of systematic errors except the uncertainty related to discretization effects. The corresponding value of the SU(2) chiral low-energy constant $\overline{\ell}_6$ is equal to $\overline{\ell}_6 = 16.2 ~ (1.0)$.


I. INTRODUCTION
The investigation of the physical properties of the pion, which is the lightest bound state in Quantum Chromodynamics (QCD), can provide crucial information on the way low-energy dynamics is governed by the quark and gluon degrees of freedom. In this respect for spacelike values of the squared four-momentum transfer, Q 2 ≡ −q 2 ≥ 0, the electromagnetic (e.m.) form factor of the pion, F π (Q 2 ), parametrizes how the pion deviates from a point particle when probed electromagnetically, thus giving insight on the distribution of its charged constituents. At momentum transfer below the scale of chiral symmetry breaking (Q 2 1 GeV 2 ) the pion form factor represents therefore an important test of non-perturbative QCD.
It is well known that for Q 2 0.5 ÷ 1 GeV 2 the experimental data on the pion form factor [1][2][3]  only. This is not too surprising in view of the fact that in the time-like region the pion form factor is dominated by the ρ-meson resonance.
An interesting issue is the quark mass dependence of the pion form factor that can be addressed by SU(2) Chiral Perturbation Theory (ChPT) known at both next-to-leading (NLO) [4] and nextto-next-to-leading (NNLO) [5] orders. The determination of F π (Q 2 ) from lattice QCD simulations provides therefore an excellent opportunity for the study of chiral logarithms. The latter are particularly important in the case of the squared pion charge radius r 2 π , i.e. the slope of pion form factor at Q 2 = 0. This means also that a controlled extrapolation to the physical point is a delicate endeavour, such that one would ideally like to perform the computation directly at the physical pion mass.
Initial studies of the pion form factor using lattice QCD in the quenched approximation date back to the late 80's [6,7] giving strong support to the vector-meson dominance hypothesis at low Q 2 . Studies of F π (Q 2 ) employing unquenched simulations have been carried out in Refs. [8][9][10][11][12][13][14][15] using pion masses above the physical one and adopting ChPT as a guide to extrapolate the lattice results down to the physical pion point. Recently a computation of F π (Q 2 ) at the physical pion mass has been provided in Ref. [16].
In this work we present a determination of the pion form factor using the gauge configurations generated in Ref. [17] by the European Twisted Mass Collaboration (ETMC) with N f = 2 twisted-mass quarks at maximal twist, which guarantees the automatic O(a)-improvement [18].
The calculations are carried out at a fixed lattice spacing (a 0.09 fm) and with pion masses equal to its physical value, 240 MeV and 340 MeV. Momentum is injected using non-periodic boundary conditions in order to get values of Q 2 between 0.01 GeV 2 and 0.25 GeV 2 . It will be shown that our data can be successfully analyzed using SU (2) ChPT at NLO without the need of the scale setting.
Our final result for the squared pion charge radius is where the error includes several sources of systematic errors except the uncertainty related to discretization effects. The corresponding value of the NLO SU(2) low-energy constant (LEC) 6 is equal to Our result (1) is obtained at a fixed value of the lattice spacing and therefore the continuum limit still needs to be evaluated. We note that discretization effects in our calculations of the pion form factor start at order O(a 2 ) (see Section III) and that our finding (1) is consistent with the experimental value r 2 exp. π = 0.452 (11) fm 2 from PDG [19]. This suggests that the impact of discretization effects on our result (1) could be small with respect to the other sources of uncertainties.
The plan of the paper is as follows. In Section II we describe the lattice setup adopted in this work, while the procedures adopted to extract the pion form factor from appropriate ratios of 3-and 2-point correlators are discussed in Section III. The lattice data for the pion form factor F π (Q 2 ) are presented in Section IV and in the Appendix, while our fitting procedures based on ChPT are described in Section V. The results of the extrapolations to the physical point and to the infinite lattice volume are collected in Section VI. Our conclusions are summarized in Section VII.

II. LATTICE ACTION
The results presented in this paper are based on the gauge configurations generated in Ref. [17] by the ETMC with Wilson clover twisted mass quark action at maximal twist [20], employing the Iwasaki gauge action [21]. The measurements are performed on N f = 2 ensembles with pion mass at its physical value, 240 MeV and 340 MeV, respectively. The lattice spacing is a 0.0914 (15) fm for all the ensembles [17]. In Table I [17].
In addition to the relevant input parameters we give the lattice volume (L/a) 3 × T /a, the number of evaluated configurations N conf and the pion mass Mπ in lattice units (with statistical error).
Both the sea and valence quarks are described by the Wilson clover twisted mass action. The Dirac operator for the light quark doublet consists of the Wilson twisted mass Dirac operator [20] combined with the clover term, namely in the so-called physical basis where D = γ µ (∇ * µ + ∇ µ )/2, ∇ µ and ∇ * µ are the forward and backward lattice covariant derivatives, and W cr = −(a/2)∇ * µ ∇ µ + m cr with m cr being the critical mass. Moreover, µ is the average up/down (twisted) quark mass, a is the lattice spacing and r = 1 the Wilson parameter. The operator D acts on a flavour doublet spinor ψ = (u, d) T . Finally, c sw is the so-called Sheikoleslami-Wohlert improvement coefficient [22] multiplying the clover term. In our case the latter is not used for O(a) improvement but serves to significantly reduce the effects of isospin breaking [17].
The critical mass has been determined as described in Refs. [23,24]. This guarantees that all physical observables can be extracted from lattice estimators that are O(a) improved by symmetry [18], which is one of the main advantages of the Wilson twisted mass formulation of lattice QCD.

III. THE PION FORM FACTOR
The pion form factor can be computed from the matrix elements of the e.m. vector current between pion states, yielding where q µ = (p µ − p µ ) is the 4-momentum transfer and Q 2 ≡ −q 2 . As detailed in Ref. [9], up to discretization effects of order O(a 2 ) it is enough to compute in Eq. (5) only the connected insertion of the single flavor currentū(x)γ µ u(x) with unitary charge.
Working in Euclidean space-time, we can access the region of space-like momentum transfer, Q 2 = −q 2 > 0, by evaluating ratios of pion two-point and three-point functions with the vector current insertion. To inject arbitrary momenta, we make use of non-periodic boundary conditions (BCs) [25][26][27] on the quark fields. Enforcing ψ(x + e i L) = e 2πiθ i ψ(x) on the quark field ψ, changes the momentum quantisation condition in finite volume to p i = 2πθ i L + 2πn i L . This is depicted in Figure 1 for the pion three-point function with independent values of the vector θ for the three quark lines. Since the ETMC gauge ensembles have been produced by imposing antiperiodic BCs in time, the same conditions are applied also to the valence quarks choosing 2πθ 0 /T = π/T . Moreover, the use of different BCs in space for sea and valence quarks produces unitarity violating finite volume effects, which are however exponentially small [28][29][30]. For the case of twisted mass quarks, this setup was first studied in Ref. [9] in the Breit frame ( p = − p), which results in a squared 4-momentum transfer independent of the pion mass, viz.
In this work the spatial components of the vector θ are always chosen to be equal each other, i.e. θ = {θ, θ, θ}.
Following Ref. [9] the required correlation functions can be evaluated efficiently through the usage of the so-called one-end-trick combined with spatial all-to-all propagators from stochastic time-slice sources and the sequential propagator method for the insertion (see Ref. [31] for the idea first applied to moments of pion parton distribution functions). Since the spatial matrix elements of the vector current are vanishing in the Breit frame, we have to compute the following correlation where is the interpolating operator annihilating the π + , t is the time distance between the vector current insertion and the source and t is the time distance between the sink and the source.
As it has been shown in Ref. [18], as well as the γ 5 -hermiticity property, one gets that: i) the correlators (6) and (7) are real, and . Thus, we have E π ( p) = E π (− p) and the discretization effects in both C 2pt (t, p) and C 3pt 0 (t, t , p, − p) start automatically at order O(a 2 ).
Taking the appropriate limits with T being the time extent of the lattice, one obtains in the Breit frame where G 2 π is the amplitude of the 2-point correlation function. Since we work from now on exclu-sively in the Breit frame, we will drop the second momentum argument and write Now, we can construct the ratio which has the following combined limit To extract F π (Q 2 ), we compute the renormalisation constant of the vector current, Z V , from the ratio of the two and three-point functions at zero momentum transfer and the known normalisation In practice, Eqs. (6-7) are evaluated by first generating stochastic sources ξ a,α r ( x, t) (r = 1, ..., N ) at a single (randomly chosen) time-slice, that for ease of notation we conventionally put in what Here, a(b) and α(β) are colour and Dirac indices, respectively, and we remind that t represents the time distance from the source. The stochastic source ξ r is manifestly zero for all t = 0. Setting one can estimate owing to γ 5 -hermiticity γ 5 D u γ 5 = D † d and γ 2 5 = 1. In Eq. (14) the pion momentum p is given by p = 2π θ/L. At fixed values of t (the time distance between the sink and the source) the so-called sequential propagator is computed as Then, one estimates the three-point function in the Breit frame kinematics from We determine F π (Q 2 ) from the ratio defined in Eq. (10) in two different ways: for the first one we compute the double ratio and we extract the pion form factor from its large time distance behavior We denote this estimate as the numerical one. The second estimate consists in replacing the pseudoscalar two-point function by its analytical expression, i.e. we fit to the data for the two-point function at large Euclidean times to determine the amplitude G π and the energy E( p). Next we define where we replace the data for the two-point function by its analytical expression using the best fit parameters. Then we calculate the double ratio from which the pion form factor can be obtained as The analytical estimate (21) may have the advantage of being less noisy than the numerical one (18), because the data for the two-point function at large t can be noisy, in particular for the largest values of | p|.
A further improvement is to replace in Eq. (19) the pion energy E π ( p), extracted from the 2-point correlator C 2pt (t, p), with the corresponding value from the dispersion relation where M π (L) is pion mass extracted from the 2-point correlator at rest. Indeed, in Figure 2 we show the measured energy levels (aE π ) 2 in lattice units as a function of the squared momentum |a p| 2 for the two ensembles cA2.30.24 and c.A2.30.48 with L/a = 24 and 48, respectively. The data is described reasonably by the dispersion relation (22) However, the use of non-periodic BCs is expected to produce further FVEs in the dispersion relation (22). Such corrections have been investigated in Ref. [32] using partially quenched ChPT at NLO, finding that the pion momentum p = 2π θ/L acquires an additive correction term 2π K/L, where the components of the vector K are given by 2 being the elliptic Jacobi function and its derivative.
For a better visualization of the effects of the additive correction (24) we consider the dimensionless quantity c 2 , defined as which in absence of FVEs on the momentum should be equal to unity. In Figure 3  In order to minimize excited state effects, the source-sink separation is fixed to t = T /2. On each gauge configuration, multiple source time slices are chosen randomly across the whole time extent, which has been shown to decorrelate measurements from different gauge configurations.
The statistical analysis is performed using the blocked bootstrap method. Before closing this Section, we address briefly the estimate of the renormalization constant of the vector current, Z V , which can be obtained form the plateau of the ratio (11). We remind that the latter one involves 2-and 3-point correlation functions with pion at rest and corresponds to fix the absolute normalization of the pion form factor, F π (Q 2 = 0) = 1. The data for the ratio (11) exhibit nice plateaux in an extended time region (t/a 5) and allow to extract Z V with a very high statistical precision (≈ 0.01%). The resulting values of Z V do depend upon the quark mass as a pure discretization effect (see also Ref. [9]). The extrapolation to the chiral limit provides therefore the value of the renormalization constant Z V , which is indeed defined in such a limit.
Using a linear fit in the (bare) quark mass 1 we get Z V = 0.6679 (1)

B. Pion Electromagnetic Form Factor
The lattice data obtained for F π (Q 2 ) as a function of Q 2 in physical units for all the gauge ensembles of Table I is shown in Figure 6 and collected in the Appendix together with the values 1 A linear dependence on the quark mass is not in contradiction with the O(a) improvement of the ratio (11), since terms proportional to a 2 µΛQCD may be dominant with respect to terms proportional to a 2 µ 2 . It is visible that the errors of our lattice data are compatible with the ones of the experimental data. In particular for small Q 2 (right panel) the errors of the lattice data are significantly smaller than the errors of the experimental single data points. Of course, the experimental points have a much denser coverage of Q 2 values. However, thanks to non-periodic boundary conditions our lattice data covers Q 2 values below the range where experimental data is available. Moreover, our lattice data at the physical pion point and at the largest volume is compatible with the experimental data within statistical errors.
We have collected in Fig. 7 the lattice data for the inverse pion form factor 1/F π versus the dimensionless variable (Q/M π ) 2 for the six gauge ensembles of Table I. It can be seen that the data for 1/F π exhibits an almost linear behavior with Q 2 , as expected from VMD arguments. Actually the solid lines in Fig. 7 represent the results of a quadratic fit in Q 2 , where we find that c π (s π ) 2 in accord with the VMD hypothesis. Moreover, for each pion mass the data for two different lattice volumes are compared in Fig. 7. It can clearly be seen that finite volume effects are relevant for M π L 3.

A. Chiral Extrapolation
Within SU(2) Chiral Perturbation Theory (ChPT) the expansion of the pion form factor F π (Q 2 ) in powers of the squared pion mass reads as where ∆F NLO π (Q 2 ) is the next-to-leading order (NLO) term and ∆F NNLO π (Q 2 ) the NNLO one.
Both are known [4,5] and the NLO term is explicitly given by where 6 is an SU(2) LEC, ξ ≡ M 2 π /(4πf π ) 2 and Let's define the slope s π and the curvature c π of the pion form factor in terms of its expansion in powers of (Q/M π ) 2 as At NLO one has which show that the LEC 6 governs only the value of the slope s NLO π . Once the value of 6 is fixed by the reproduction of the experimental value of the pion charge radius, i.e. 6 14.6 (see Ref. [9]), it turns out that c NLO π (s NLO π ) 2 , which is in contradiction with the VMD phenomenology observed both in the experimental data and in our lattice results up to a pion mass of 340 MeV (see Fig. 7). In Ref. [9] it was found that, using 6 14.6 the NLO term (28) works only for very low values of both Q 2 (Q 2 0.03 GeV 2 ) and the pion mass (M π 300 MeV). Effects from NNLO and higher order terms in the chiral expansion (27) become more and more important as the value of Q 2 increases. In particular, the curvature c π is found to be almost totally dominated by NNLO effects [9]. The latter however depend on several LECs (see Ref. [5]).
In order to avoid the need of many LECs let's consider the inverse of the pion form factor.
Using Eq. (27) the SU(2) ChPT expansion of 1/F π (Q 2 ) reads as where on the r.h.s. the term in the square brackets represent the NNLO correction. Because of the observed VMD phenomenology (see Fig. 7), the NNLO term ∆F NNLO π (Q 2 ) in Eq. (33) is expected to be almost compensated by the square of the NLO one ∆F NLO π (Q 2 ), leading to a small residual NNLO correction in the inverse pion form factor. This means that 1/F π (Q 2 ) is dominated by the NLO approximation at least in the range of values of Q 2 and M π covered by our simulations, i.e. Q 2 0.25 GeV 2 and M π 340 MeV. Thus, we can profit from the above feature by using the following ansatz for the chiral extrapolation of the inverse pion form factor where the last term in the r.h.s. parametrizes NNLO effects, which we stress are expected to be small. Eq. (34) depends only on three unknowns, namely 6 , b 1 and b 2 , which we determine by fitting our data.

B. Finite Volume Effects
As illustrated in Fig. 7, our data for the pion form factor suffer from finite volume effects (FVEs). In this work we follow three strategies to correct for FVEs, profiting from the two lattice volumes available at each value of the quark mass.
The three strategies are as follows: A) make use of the SU(2) ChPT prediction derived at NLO in the Breit frame [33,34]. The correction factor K FVE (Q 2 , L) reads explicitly where C is a parameter to be determined in the fitting procedure, Q 2 = 4(2π θ/L) 2 and with Θ(τ, θ) ≡ ∞ n=−∞ e −τ (n+θ) 2 being the elliptic Jacobi function.
B) use a phenomenological ansatz, inspired by the asymptotic expansion of Eq. (37), given by where C 1 and C 2 are parameters to be determined in the fitting procedure.
C) use only the largest volume available at each pion mass and assume that FVEs are negligible for these volumes (i.e. putting K FVE = 0).
We want to point out that our fitting ansatz (36) is defined in terms of dimensionless quantities only, namely ξ, Q 2 /M 2 π and M π L, and therefore the knowledge of the lattice scale is not required.

VI. EXTRAPOLATIONS TO THE PHYSICAL POINT
In this Section we perform the chiral and infinite volume extrapolations of the lattice data adopting our fitting ansatz (36). Various sources of systematic effects have been taken into account, namely • the numerical and analytical estimates of the pion form factor given by Eqs. (18) and (21), respectively. The corresponding uncertainty will be denoted by () ratio ; • the time extensions [t min , T /2 − t min ] chosen for the plateaux of the double ratios (17) and (20) corresponding to t min = 10 and 12. The corresponding uncertainty will be denoted by () f it−range ; • either the inclusion of all the six gauge ensembles of Table I or the restriction to the two gauge ensembles cA2.09.XX at the physical pion mass. The corresponding uncertainty will be denoted by () Mπ ; • either the inclusion (b 1 = 0 and b 2 = 0) or the exclusion (b 1 = b 2 = 0) of the NNLO effects in Eq. (34). The corresponding uncertainty will be denoted by () ChP T ; • the FVEs evaluated according to the three procedures A, B and C, described in Section V B.
The corresponding uncertainty will be denoted by () F V E ; • the inclusion of all Q 2 values or the restriction to Q 2 ≤ 2M 2 π . The corresponding uncertainty will be denoted by () Q 2 −range .
The quality of our fitting procedure is illustrated in Fig. 8  KFVE(Q 2 , L) given by Eqs. (34) and (37), respectively, for all the six ensembles of Table I. The bands correspond to statistical uncertainties only.
infinite volume limit, are compared with the experimental data in Fig. 9. Eq. (28) of Ref. [35] As far as the pion charge radius is concerned, our fitting ansatz (36) implies that at the physical pion point and the infinite volume limit one has Thus, our final result at a fixed lattice spacing (a 0.09 fm) reads which is consistent with the experimental value r 2 exp. π = 0.452 (11) fm 2 [19]. This suggests that the impact of discretization effects on our result (1) could be small with respect to the other sources of uncertainties.
The value of the NLO SU(2) LEC 6 , appearing in Eq. (40) and corresponding to our result (41), is equal to

VII. SUMMARY AND DISCUSSION
We have presented an investigation of the electromagnetic pion form factor, F π (Q 2 ), at small values of the four-momentum transfer Q 2 ( 0.25 GeV 2 ), based on the gauge configurations generated by ETMC with N f = 2 twisted-mass quarks at maximal twist including a clover term.
Momentum is injected using non-periodic boundary conditions and the calculations are carried out at a fixed lattice spacing (a 0.09 fm) and with pion masses equal to its physical value, 240 MeV and 340 MeV. We have successfully analyzed our data using Chiral Perturbation Theory at nextto-leading order in the light-quark mass. For each pion mass two different lattice volumes are used to take care of finite size effects. Our final result for the squared charge radius is r 2 π = 0.443 (29) fm 2 , where the error includes several sources of systematic errors except the uncertainty related to discretization effects. The corresponding value of the SU(2) low-energy constant 6 is equal to 6 = 16.2 (1.0). Our result is consistent with the experimental value r 2 exp. π = 0.452 (11) fm 2 [19] as well as with other lattice estimates (see Ref. [36]). This suggests that the impact of discretization effects on our result could be small with respect to the other sources of uncertainties.  [37][38][39], Lemon [40], DDαAMG [41] and R [42] have been used.

APPENDIX
In Tables II-VII we collect the values adopted for the vector θ = {θ, θ, θ}, the squared 3momentum | p| 2 = 4π 2 | θ| 2 /L 2 in lattice units, the squared 4-momentum transfer Q 2 = −q 2 = 4| p| 2 in units of the pion mass and the values of the pion form factor F π (Q 2 ) in the case of the six ensembles of Table I.