On a novel evalutation of the hadronic contribution to the muon's $g-2$ from QCD

We evaluate the hadronic contribution to the $g-2$ of the muon by deriving the low-energy limit of quantum chromodynamics (QCD) and computing in this way the hadronic vacuum polarization. The low-energy limit is a non-local Nambu--Jona-Lasinio (NJL) model that has all the parameters fixed from QCD, and the only experimental input used is the confinement scale that is known from measurements of hadronic physics. Our estimations provide a novel analytical alternative to the current lattice computations and we find that our result is close to the similar computation performed from experimental data. We also comment on how this analytical approach technique, in general, may provide prospective estimates for hadronic computations from dark sectors and its implication in BSM model-building in future.


I. INTRODUCTION
Since the original computation for the electron from first principles [1], originating from Dirac equation, the lepton anomalous magnetic moments continue to be very important observables for precision tests of the Standard Model (SM) [2]. Recent data seems to indicate a tension with the theoretical prediction for the anomalous magnetic moment of the muon with the recent experimental value for the anomalous magnetic moment of the muon is [3,4] g µ /2 = 1 + a µ = 1.001 165 920 8 (6). (1) The Particle Data Group (PDG) gives an updated value for the muon anomaly in the form [5] a exp µ = 116 592 091(54)(33) × 10 −11 .
This precision clearly is a challenge for the theoretical side to increase the precision of the prediction [4] The theoretical results for the muon anomalous magnetic moment in the SM are traditionally represented as a sum of three parts, with a QED µ , a EW µ being the leptonic and electroweak parts, respectively, and a had µ is the contribution involving the electromagnetic currents of quarks.
The hadronic part a had µ in the SM is related to quark contributions to the electromagnetic currents.
The difference ∆a µ = a exp µ − a SM µ = 268(63)(43) × 10 −11 (7) could be due to new unknown physics beyond the SM, but it is not statistically significant off yet [6], the main idea behind this being that contributions from unknown virtual particles not part of the SM might enter the calculations.
In general, theoretical estimates are very precise for what one should expect from quantum electrodynamics (QED), but fall short in the case of the hadronic contributions, due to the known difficulties to treat quantum chromodynamics (QCD) at low energies. The general accepted technique is to use experimental results from e + e − scattering into hadrons measured in collider experiments [6]. Two key ingredients of this contribution are the hadronic vacuum polarization (HVP) and the high-order hadronic light-by-light scattering (HLbL). Of these two contributions, the former is the most critical one, due to the current inability to compute this contribution starting right from the Lagrangian of QCD. Some recent evaluation of the HVP from experimental data is given in Ref. [7][8][9] [10], showing that the HVP correction they obtain moves the ballpark of the muon g − 2 value back into the SM field. In turn, this would imply that the technique using experimental values from the colliders probably underestimates this contribution.
Working on QCD, one generally makes the use of effective models. However, it is often unknown if such models could be straightforwardly obtained from the original Lagrangian.
Still, successful results have been obtained from some of these models. In the very early days of the study of the muon g − 2 problem, attempts were made to derive the HVP contribution from such effective models like for instance the Nambu-Jona-Lasinio (NJL) model [11], detailed by de Rafael in Ref. [12] and more recently [13]. However, due to the large set of undetermined parameters entering in such effective theories, this kind of approach in this early, primitive stage was abandoned in favor of the use of experimental data and lattice QCD calculations.
Inspired by such an approach, in this article we will show how an effective field theory can be derived from QCD, starting directly from the Lagrangian level. The model is a non-local Nambu-Jona-Lasinio model, having all the parameters properly fixed. A first attempt in this direction was given in Ref. [14] in order to determine the proper low-energy limit of the theory 1 . In this work, we fix an error in this publication and show how the effective NJL model comes out naturally from QCD. Based on these first principles, we will evaluate the HVP contribution to the muon (g − 2).

II. BASIC EQUATIONS FOR NJL MODEL
Our starting point is the well-known QCD Lagrangian with covariant derivative D µ = ∂ µ + igT a A a µ and the field strength tensor components F a µν defined by igT a F a µν = [D µ , D ν ]. The sum over i is understood to run over the quark flavors and colors. Throughout this paper we work with the Minkowskian metric g µν = diag(1, −1, −1, −1). Calculating the Euler-Lagrange equations, one obtains These classical equations of motion are the starting point for a tower of Dyson-Schwinger equations. In order to study these equations, we use the method proposed by Bender, Milton and Savage [15], details of which can found in Refs. [14,[16][17][18][19]. For the purpose of this publication we sketch the main steps here, skipping contributions from BRST ghosts for simplicity.
Enlarging the Lagrangian of the classical action by adding corresponding source terms A a µ J µ a ,q i η i andη i q i , one obtains the exponential of the generating functional. Functional 1 Unfortunately, this publication contained a mistake that made the conclusions unreliable.
derivatives of this generating functional lead to the Dyson-Schwinger analogue of the Euler-Lagrange equations, expressed in terms of Green functions for the fields. The set of equations in Landau gauge ξ = 0 we start with is given by where the one-, two-and three-point Green functions are given by A 1a

y) . The expected solutions can be written in the form
where η a µ are the coefficients of the polarization vector with η a µ η µ b = δ ab , φ(x) is a scalar field and ∆(x − y) is the propagator of the scalar field. The three-point function can be set to zero. For the one-point functions we obtain Using , the first differential equation (12) takes the form In the 't Hooft limit N c → ∞, λ := N c g 2 ≫ 1 finite but large, this set of equations yields a Nambu-Jona-Lasinio model in a straightforward way. Indeed, for this case, we can perform a perturbation series expansion φ(x) = φ 0 (x) + φ 1 (x) + O(g 2 ) in g, obtaining at leading order while the next-to-leading order yields A. Zeroth order solution and Green function Note that ∆(0) is a constant. Therefore, m 2 = 2λ∆(0) can be considered as the mass square of the scalar field. The leading order differential equation λφ 3 0 (x) = 0 is nonlinear, but a solution in terms of Jacobi's elliptic functions exists, with where µ and θ are integration constants. sn(z|κ) is Jacobi's elliptic function of the first kind.
Given this solution, the second differential equation can be solved by noting that is solved by a Green function written in momentum space as with andẐ The mass spectrum is given by At this point the circle for the mass of the scalar field is closed. Inserting back the Fourier transform of the propagator (19) into m 2 = 2λ∆(0) results in This self-consistency equation provides the proper spectrum of a Yang-Mills theory with no fermions [19], in very close agreement with lattice data.

B. First order solution
The convolution of the propagator ∆ with the right hand side of Eq. (15) leads to The first term renormalizes the fermion mass and can taken to be zero by choosing the renormalization condition S(0) = 0. The second term yields a Nambu-Jona-Lasinio (NJL) interaction in the equation of motion of the quark.
Inserting φ(x) into the Dirac equation (12), in the 't Hooft limit the term φ 0 is negligible small compared to the NJL interaction term φ 1 . This can be realized by noting that φ 0 ∼ λ 1/4 while φ 1 ∼ λ. In the strong coupling limit λ ≫ 1, for the quark one-point function we have to retain only the NJL term. Therefore, we obtain The Note that η η a µ η b ν = δ ab g µν , where η symbolizes the polarizations. In addition, one traces out the color degrees of freedom with tr(T a T a ) = N c C F , C F = (N 2 c − 1)/(2N c ), and are spinors in Dirac and flavor space, only. This leads us to the NJL lagrangian The Fierz rearrangement of the quark fields yields

C. Bosonization
Let Γ α be a set of Dirac and flavor matrices containing not only the Dirac structures 1, iγ 5 , γ µ and γ µ γ 5 from the Fierz rearrangement but also the flavor matrices 1l and 1 2 λ α relating quarks of equal and different flavor i and j in adjoint representation. Γ α obeys the conjugation rule γ 0 Γ † α γ 0 = Γ α , where α denotes the components of the adjoint flavor representation. Accordingly, the spinor ψ(x) spans over all these spaces. The most prominent degrees of freedom are the scalar-isoscalar and pseudoscalar-isovector degrees which can formally be combined as four vector. As the coefficients of these two contributions are the same, one can reinterpret the sum over these 1 + 3 = 4 degrees of freedom as a sum over four-vector components. The next step is to apply the bosonization procedure exemplified in Ref. [22] by adding scalar-isoscalar and pseudoscalar-isovector mesonic fields as auxiliary fields M α (w) = (σ(w); π(w)) at an intermediate space-time location w = (x + y)/2, coupled to the nonlocal fermionic currents. The result of the Fierz rearrangement can be expressed as NJL action (G = 2 d 4 z∆(z)). By performing a nonlocal functional shift the nonlocal quartic fermionic interaction can be removed. Instead, the fermion field starts to interact nonlocally with the mesonic fields, After Fourier transform, in momentum space one obtains where the symbols with tilde are used for the Fourier transformed quantities. The final step in the bosonization is to integrate out the fermionic fields, in the general case leading to [22] S bos = − N c g 2 4G where det denotes the direct product of a functional and an analytical determinant, the former in the Fock space transition between space-time points x and y, the latter in the Dirac and flavor indices.

D. Mean field approximation
Expanding the bosonic fields σ(x) =σ+δσ(x) and π(x) = δ π(x) about the vacuum expectation valueσ = σ , the zeroth order expansion coefficient is the mean field approximation, leading to the simplified NJL action After Fourier transform, in momentum space one has with the unit space-time volume V (4) , where (G = 2∆(0)) is the dynamical mass of the quark. The bosonization leads to On the other hand, one has ln det(p / − M q (p)) = tr ln(p / − M q (p)) = 1 2 4N f ln p 2 − M 2 q (p) . The quantityσ can be determined by variation of the action S bos with respect to this quantity. Taking into account the dependence of M q (p) onσ, one obtains Finally, this result can be re-inserted to Eq. (36) to obtain the dynamical mass equation A similar gap equation for the g − 2 problem was shown in Ref. [13]. In this article, however, we derived the gap equation directly from the QCD Lagrangian.

III. SOLVING THE GAP EQUATION
At this point we can insert∆(p) from Eq. (19) into Eq. (39) in order to obtain the gap equation for the dynamical quark mass -or to be more precise the couple of gap equations, if taking into accoung Eq. (23) as well. However, in order to make the calculation feasible, we recognize that the dependence on the mass m of the scalar field is subdominant, and this mass can be neglected compared to the mass of the quark. For m = 0 one has κ = −1, (40) m n = (2n + 1) 2p 2 /2K(−1) = (2n + 1)m 0 is the glue ball spectrum, with the ground state given by m 0 = m G (−1) = 2p 2 /2K(−1) and K(z) is the complete elliptic integral of the first kind. As a further simplification we calculate the dynamical quark mass at zero momentum, p = 0. In this case we obtain (41) where we have performed a Wick rotation to the Euclidean domain. As this integral is UV singular, we integrate the momentum up to a cut Λ to obtain where we have used the dimensionless quantities x = m 0 /Λ and y = M q /Λ, assuming that M q ≪ Λ. Reinserting into Eq. (41) leads to the gap equation where κ = N f N c /π and α s = g 2 /4π. We note that the cut-off completely disappeared except for the ratio m q /Λ that, for the light quarks, is negligible small.
For the QCD cut-off Λ = 1 GeV, the average mass of the u and d quarks is taken to be m q = 0.003415(48) GeV [5]. The ground state of the glue ball spectrum is given by the f 0 (500) resonance, measured as m 0 = 0.512(15) GeV [23]. Using N c = 3, N f = 6 and α s (3.1 GeV) = 0.256506 we obtain M q = 0.427(29) GeV.

IV. HADRONIC VACUUM POLARIZATION
Inspired by the approach in Ref. [12], next we will evaluate the contribution to the hadronic vacuum polarization, assuming that a NJL approximation holds [14,17]. Looking at the Fierz decomposition as shown in Eqn. (28), one obtains where which agrees well with the analysis in the preceding section, provided we evaluate the gap equation as in Eq. (43). Using Ref. [12], we evaluate The coefficient P 1 determines the contribution called "had 1a". It is defined by where Π (H) and Π (1) Γ(n, ε) is the incomplete gamma function, but Γ(1, ε) is an analytic expression, Using these formulas, we obtain a u,d µ (had 1a) = 452(67) · 10 −10 .
In order to have a clearer understanding of the meaning of this result, we present also the strange quark contribution. This will yield a s µ (had 1a) = 232(34) · 10 −10 .
The overall is The error bar is not yet competitive to decide if BSM physics is needed but nevertheless in closed agreement with the experimental value as obtained in [6][7][8] from experiments in hadron physics.
Finally, we want to analyze the contribution to the error due to the choice of the 't Hooft limit: Ng 2 constant and N → ∞. There have been several studies on lattice to estimate the error of such an approximation ( [24,25] and references therein). The main conclusion is that the next-to-leading order correction to any observable goes like being c 1 = O(1), a numerical factor. This same pattern is seen in the spectrum of a Yang-Mills theory without quarks where, for the ground state, one sees [26] m 0 ++ √ σ = 3.28(8) + 2.1(1.1) where σ is a mass scale proper to strong interactions and obtained by experiment. So, this can be estimated of the same magnitude as the error we obtained from QCD data at worst.

V. CONCLUSIONS AND OUTLOOK
To summarise, using technique devised by Bender, Milton and Savage, in Ref [15] the Dyson-Schwinger equations for quantum chromodynamics in differential form was revisited. Following Ref. [15], in this article we discussed the hadronic contributions to the muon anomalous magnetic moment following NJL model as the low energy effective theory description of QCD, as shown in Eq. (26). We provided a full derivation of the HVP contribution to the anomalous magnetic moment a = (g − 2)/2 of the muon from first principles, starting from the QCD partition function and the effective mass for the quarks as shown in Eq. (39).
Our result as obtained in Eq. (52) is in close agreement with the Muon g − 2 Theory Initiative [6], as obtained from experimental data in Ref. [7][8][9]. In doing so, we have shown a possible new analytical approach as an alternative to lattice calculations. Our approach provides a theoretical framework for the application of QCD to several other applications and the opportunity to investigate future studies model-building for BSM physics in the dark sector just by using analytical methods. The next step will be to include other quark flavors which is beyond the scope of the current manuscript. Moreover, following the same approach and using NJL model as the low energy EFT for QCD, we also can perform a complete proof of confinement in QCD in our future studies.
We hope to improve our computations in the near future to reduce the error bar significantly.