{'fraction_non_alphanumeric': 0.04161563221010055, 'fraction_numerical': 0.01449853305608568, 'mean_word_length': 4.124507371186688, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 16, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Vector and axial vector mesons in a nonlocal chiral quark model

Basic features of nonstrange vector and axial vector mesons are analyzed in the framework of a chiral quark model that includes nonlocal four fermion couplings. Unknown model parameters are determined from some input values of masses and decay constants, while nonlocal form factors are taken from a fit to lattice QCD results for effective quark propagators. Numerical results show a good agreement with the observed meson phenomenology.


I. INTRODUCTION
Given the nonperturbative character of quantum chromodynamics (QCD) in the low-energy regime, the analysis of hadron phenomenology starting from first principles is still a challenge for theoretical physics. Although substantial progress has been achieved in this sense through lattice QCD (LQCD) calculations, this approach shows significant difficulties, e.g. when dealing with small quark masses or with hadronic systems at nonzero chemical potentials. Thus it is important to study the consistency between the results obtained through lattice calculations and those arising from effective models for strongly interacting particles. For two light flavors it is believed that QCD supports an approximate SU(2) chiral symmetry that is dynamically broken at low energies, where pions play the role of the corresponding Goldstone bosons. The well-known Nambu−Jona-Lasinio (NJL) model [1,2], in which light mesons are described as fermion-antifermion composite states, is a simple effective approach that shows these features. In the NJL model quarks interact through a local four-fermion coupling, leading to relatively simple Schwinger-Dyson and Bethe-Salpeter equations. Now, as a step toward a more realistic approach to low-energy QCD, it is worth it to consider extensions of the NJL model that include nonlocal interactions [3]. In particular, this is supported by lattice calculations, which lead to a given momentum dependence of both the mass and the wave function renormalization (WFR) in the effective quark propagators [4,5]. It is also seen that nonlocal extensions of the NJL model do not exhibit some problems that are present in the local theory. For example, nonlocal interactions regularize the model in such a way that the effective interaction is finite to all orders in the loop expansion, thus model predictions are less dependent on the parameterizations, and there is no need to introduce extra cutoffs [6].
Previous works on nonlocal NJL-like (nlNJL) models, focused on different aspects of strong interaction physics, can be found in the literature. These include the study of vacuum hadronic properties considering either two [7][8][9][10][11][12][13][14] or three [15] active quark flavors, and various nonlocal form factor shapes. In addition, this framework has been used to describe the chiral restoration transition for hadronic systems at finite temperature and/or chemical potential (see e.g. Refs. [16][17][18][19][20][21]). In this work, following the proposal in Refs. [11,13], we consider a model in which nonlocal form factors lead to a momentum dependence of the mass and WFR in the quark propagator, hence the actual shape of these form factors can be taken from the data obtained through lattice calculations [13,19]. We concentrate here in particular in the incorporation of explicit vector and axial vector interactions. Therefore, besides the previously considered couplings between scalar and pseudoscalar quark-antiquark currents, in our model we include couplings between vector and axial vector nonlocal currents satisfying proper QCD symmetry requirements. In fact, nonlocal models including vector and axial vector currents have been previously considered in Ref. [9].
However, those models do not include a momentum-dependent WFR of quark propagators, which is required in order to perform the comparison with lattice QCD results. We dedicate the first part of the paper to work out the formalism in order to derive analytical expressions for some basic vector meson properties, such as masses and decay parameters. Then we present numerical results obtained by taking the nonlocal form factors from a fit to lattice QCD data. It is seen that, after fixing unknown coupling constants so as to reproduce some input meson observables, the model provides an adequate phenomenological description of the considered vector meson properties.
The article is organized as follows. In Sect. 2 we introduce the model and derive the corresponding gap equations at the mean field level. In Sect. 3 we describe the vector meson sector, obtaining analytical results for meson masses and decay amplitudes. The numerical and phenomenological analyses are included in Sect. 4, while in Sect. 5 we present a summary of our work. Finally, in Appendixes A and B we collect some analytical expressions and describe the calculation procedure.

II. MODEL
We consider a two-flavor chiral quark model that includes nonlocal vector and axial vector quark-antiquark currents. Since our aim is to choose form factors that are in agreement with LQCD calculations, it is convenient to work in Euclidean space, where nonlocal interactions are well defined [3]. The corresponding effective action is given by is the current quark mass matrix. We will work in the isospin symmetry limit, assuming m u = m d , which will be called from now on m c . The fermion currents are given by [13] where τ a , a = 1, 2, 3, are the Pauli matrices, while u( Eqs. (2) include the usual scalar (I = 0) and pseudoscalar (I = 1) quark-antiquark currents [11,12], as well as vector and axial-vector quark-antiquark currents that transform as either isospin singlets or triplets. In addition, we consider a coupling between "momentum" currents j M (x) [11,13], which involve derivatives of the fermion fields. The presence of this interaction is naturally expected as a correction arising from the underlying QCD dynamics. Whereas in a local theory, at the mean field level, it would simply lead to a redefinition of fermion fields, in our nonlocal scheme it leads to a momentum-dependent wave function renormalization of the quark propagator, in consistency with LQCD analyses. For convenience, we have chosen to take a common coupling constant G S for both the scalar/pseudoscalar and momentum quark interaction terms. Notice, however, that the relative strength between these terms is controlled by the mass parameter κ in j M (x). Finally, the functions f (z), g(z), h(z), h 0 (z) and h 5 (z) are covariant form factors responsible for the nonlocal character of the interactions. Notice that, in order to guarantee chiral invariance, the form factor g(z) has to be equal for the scalar and pseudoscalar currents j S (x) and j a P (x), and the same applies to the form factor h(z) entering the vector and axial vector currents j a V µ (x) and j a Aµ (x). To work with mesonic degrees of freedom, we proceed to perform a bosonization of the fermionic theory [3]. This is done in a standard way by considering the corresponding partition function Z = DψDψ exp[−S E ], and introducing auxiliary bosonic fields σ 1 (x), σ 2 (x) [scalar, related respectively to the currents j S (x) and j M (x)], π a (x) (pseudoscalar), v 0 µ (x), v a µ (x) (vector) and a 0 µ (x), a a µ (x) (axial vector), where indices a run from 1 to 3. After integrating out the fermion fields the partition function can be written as where S bos E stands for the Euclidean bosonized action. In momentum space, the latter is given by where the operator A(p, p ′ ) reads where and the trace extends over Dirac, color, flavor and momentum spaces. In the same way, the bosonized effective action in Eq. (4) can be expanded in powers of meson fluctuations as where the mean field action per unit volume reads [13] S MFA the trace acting just over Dirac space. From Eq. (7), the mean field effective quark propagator D 0 (p) is given by where the functions m(p) and z(p) -momentum-dependent effective mass and WFR-are related to the nonlocal form factors and the vacuum expectation values of the scalar fields by The mean field valuesσ 1,2 can be found by minimizing the mean field Euclidean action. This leads to the set of coupled gap equations [13] where we have defined D(p) = p 2 + m(p) 2 . The chiral quark condensates -order parameters of the chiral restoration transition-are given by the vacuum expectation values qq , where q = u, d.
The corresponding expressions can be obtained by differentiating the MFA partition function with respect to the current quark masses. Away from the chiral limit, this leads in general to divergent integrals. Since one is interested in the description of the nontrivial vacuum properties arising from strong interactions, it is usual to regularize these integrals by subtracting the free quark contributions (see e.g. Refs. [9,11,17,18]). One gets in this way

III. MESON MASSES AND DECAY CONSTANTS
We are interested in the description of vector meson phenomenology, which requires going beyond the MFA. In this section we derive analytical expressions to be used for the calculation of basic measurable phenomenological quantities, such as meson masses and decay constants. It is important to notice that pion observables, already calculated within this framework in previous works [11,13,22], need to be revisited owing to the mixing between π and a µ fields.

A. Meson masses and mixing
In general, meson masses can be obtained from the terms in the Euclidean action that are quadratic in the bosonic fields. When expanding the bosonized action we obtain where the functions G M (p 2 ), M = σ, σ ′ , π, . . . are given by one-loop integrals arising from the fermionic determinant in the bosonized action. In the case of the σ 1 , σ 2 sector the expression in Eq. (14) is given in terms of the fields σ and σ ′ , which are defined as linear combinations of σ 1 and The mixing angles θ and θ ′ are fixed in such a way that there is no σ − σ ′ mixing terms at the level of the quadratic action for p 2 = −m 2 σ (′) , where the minus sign is due to the fact that the action is given in Euclidean space. Once cross terms have been eliminated, the functions G M (p 2 ) stand for the inverses of the effective meson propagators, thus scalar meson masses are obtained by solving the equations G σ (′) (−m 2 σ (′) ) = 0. Explicit expressions for the functions G σ (′) (p 2 ) can be found in Ref. [13].
To analyze the vector meson sector one has to take into account the tensors G µν v , G µν a , G µν 0 and G µν 5 . From the expansion of the fermionic determinant we obtain where with q ± = q ± p/2. The functions G ρ,a 1 (p 2 ) and L ± (p 2 ) correspond to the transverse and longitudinal projections of the vector and axial vector fields, describing meson states with spin 1 and 0, respectively. Thus the masses of the physical ρ 0 and ρ ± vector mesons (which are degenerate in the isospin limit) can be obtained by solving the equation In addition, in order to obtain the physical states, the vector meson fields have to be normalized where Here g ρqq can be viewed as an effective ρ meson-quark effective coupling constant. Regarding the isospin zero channels, it is easy to see that the expressions for G µν 0 (p 2 ) can be obtained from those for G µν v (p 2 ), just replacing G V → G 0 and h(q) → h 0 (q). In this way, one can define for the ω vector meson a function G ω (p 2 ), obtaining the ω mass and wave function renormalization as in Eqs. (19) and (21). Similar relations apply to the axial vector sector, where G µν 5 (p 2 ) can be obtained from G µν a (p 2 ) by replacing G V → G 5 and h(q) → h 5 (q). The lightest physical state associated to this sector (quantum numbers I = 0, J P = 1 + ) is the f 1 axial vector meson, hence we denote by G f 1 (p 2 ) the form factor corresponding to the transverse part of G µν 5 (p 2 ). In the case of the pseudoscalar sector, from Eq. (14) it is seen that there is a mixing between the pion fields and the longitudinal part of the axial vector fields [23,24]. The mixing term includes a loop function G πa (p 2 ), while the term quadratic in δπ is proportional to the loop function G π (p 2 ). These functions are given by where once again we have used the definitions q ± = q ± p/2. The physical states˜ a µ and˜ π can be now obtained through the relations [23,24] where the mixing function λ(p 2 ), defined in such a way that the cross terms in the quadratic expansion vanish, is given by The pion mass can be then calculated from Gπ(−m 2 π ) = 0, where while the pion WFR can be obtained from In the case of the a 1 axial vector mesons (I = 1 triplet), since the transverse parts of the a b µ fields do not mix with the pions, the corresponding mass and WFR can be calculated using relations analogous to those quoted for the vector meson sector, namely Eqs. (19) and (21), with G a 1 (p 2 ) given by Eq. (17).

B. Pion weak decay
By definition the pion weak decay constant f π is given by the matrix elements of axial currents between the vacuum and the physical one-pion states, evaluated at the pion pole. To determine the axial currents, we "gauge" the effective action S E , introducing external gauge fields. In general, for a local theory, this is carried out just by replacing where G µ is the corresponding gauge field. In our model, due to the nonlocality of the interactions, the gauging procedure requires the introduction of gauge fields not only through the covariant derivative in Eq. (28) but also through a parallel transport of the fermion fields in the nonlocal currents (see e.g. Refs. [3,8,12]): Here x and z are the variables in the definitions of the nonlocal currents in Eq. (2), while the where s runs over an arbitrary path connecting x with y. In the case of the axial current we introduce the axial gauge fields W a µ (x), taking In addition, notice that if the action is written in terms of the original states π b and a b µ , in order to calculate the matrix element in Eq. (27) one has to take into account the mixing described in the previous subsection. Once the gauged effective action is built, the matrix elements can be obtained by taking derivatives with respect to the gauge and the physical pion fields, π, a 1 π, a 1 The resulting one-loop contributions are diagrammatically schematized in Fig. 1. Tadpole-like diagrams, which are not present in the local NJL model, arise from the occurrence of gauge fields in Eqs. (29). We finally obtain where It is important to notice that the result for f π does not depend on the path chosen for the transport function in Eq. (30) [see the comment after Eq. (42) below]. In the absence of vector meson fields, the mixing term in Eq. (33) vanishes and our expression reduces to that previously quoted in Ref. [13].
C. ρ meson-photon vertex and ρ electromagnetic decay constant Another important quantity to be studied is the ρ-photon vertex. In our nonlocal model, mesonphoton couplings receive in general contributions from the parallel transport in Eq. (29), therefore we find it important to check that the conservation of the vector current is satisfied. In addition, from this vertex we can obtain a prediction for the electromagnetic ρ → e + e − decay amplitude.
The ρ-photon vertex is given by the matrix element of the electromagnetic current between a vector meson state and the vacuum, To calculate this matrix element one can follow the procedure discussed in the previous subsection, taking now where e is the proton charge and Q = diag(2/3 , −1/3).
Once again it is possible to distinguish two contributions to Π a µν , namely Π arising from a two-vertex and a tadpole-like diagram, respectively (see Fig. 2). We obtain Here we have defined, for a given function f (p), with where s runs over a path connecting the origin with a point located at z.
which leads to It is also worth noticing that the integral in Eq.
A similar cancellation has been found in Ref. [9], where a nlNJL model that includes vector mesons without quark WFR is considered.
Let us now concentrate on the ρ electromagnetic decay constant f v , which can be defined from ρ 0 → e + e − decay: where α = e 2 /(4π) is the electromagnetic fine structure constant. It can be seen that f v is related to the trace of Π 3 µν (p) through To evaluate the transverse part of the tensor Π 3 µν we take a straight line path for the integral over s µ in Eq. (30). This leads to where f ′ (p) denotes the derivative of f with respect to p 2 . After some algebra, we obtain where respectively, while the functions α + f (q, p) have been defined as D. π 0 → γγ decay Let us analyze in the context of our model the anomalous decay π 0 → γγ. As it is well known, in the NJL model this decay is problematic: in order to reproduce the experimentally observed result it is necessary to perform quark loop momentum integrations up to infinity instead of following the cutoff prescription of the model [25]. In our framework, taking into account the discussion of gauge interactions in the previous subsections, the decay amplitude can be calculated from the matrix element In principle there are several diagrams that contribute to the amplitude at the level of one loop.
As in the case of the pion decay constant f π , since the physical π 0 stateπ 3 (p) is a combination of π and a µ fields, one has to consider the linear expansion of the bosonized action in π and in a µ . The diagrams leading to nonzero contributions are those depicted in Fig. 3. If the outgoing photons are assumed to be in states of four-momenta k 1 and k 2 with polarization vectors ε and ε (λ 2 ) ν (k 2 ), respectively, the decay amplitude can be written as where the form factorF (k 1 , k 2 ) is given by the sum of π and a µ contributions to theπ 3 state, with p = k 1 + k 2 .
The first term in the brackets, corresponding to the diagram in Fig. 3a, has been calculated (apart from an isospin factor) in Ref. [22]. One has where A(q, k 1 , k 2 ) = 1 z(q) On the other hand, the form factor F a (k 1 , k 2 ) arises from the sum of the contributions correspond-  ing to the diagrams in Figs. 3b and 3c. Although these turn out to be separately divergent, it is seen that divergent pieces cancel out and the sum is finite. We obtain where B(q, r, s) = 1 z(q) Finally, after phase space integration and sum over outgoing photon polarizations, the π 0 → γγ decay amplitude is given by Since photons are on-shell, from Lorentz invariance it is seen thatF (k 1 , k 2 ) can only be function of the scalar product (k 1 · k 2 ) = −m 2 π /2.

E. ρ → ππ decay
In general, various transition amplitudes can be calculated by expanding the bosonized action to higher orders in meson fluctuations. In this subsection we concentrate in the processes ρ 0 → π + π − and ρ ± → π ± π 0 , which are responsible for more than 99% of ρ meson decays. The decay amplitudes M(v a (p) → π b (q 1 )π c (q 2 )) are obtained by calculating the corresponding functional derivatives of the effective action, which can be written in terms of two form factorsF ρππ (p 2 , q 2 1 , q 2 2 ) and G ρππ (p 2 , q 2 1 , q 2 2 ): Only the transverse piece, driven by the form factorG ρππ (p 2 , q 2 1 , q 2 2 ), contributes to ρ → ππ decay widths. Indeed, in the isospin limit, one has where g ρππ ≡G ρππ (−m 2 ρ , −m 2 π , −m 2 π ). The form factorG ρππ (p 2 , q 2 1 , q 2 2 ) arises from the effective vertexρππ, whereρ andπ are renormalized states. Since we expand the effective action in Eq. (4) in powers of the unrenormalized fields, it is convenient to write the effective vertex in terms of the original fields ρ, π and a µ [the latter has to be taken into account due to the π − a mixing given by Eq. (23), as mentioned in previous subsections]. In this way, the form factor receives contributions from the diagrams sketched in Fig. 4. One has where G ρππ (p 2 , q 2 1 , q 2 2 ), G ρπa (p 2 , q 2 1 , q 2 2 ) and G ρaa (p 2 , q 2 1 , q 2 2 ) are one-loop functions that arise from the expansion of the effective action. The explicit forms of these functions, which can be obtained after a rather lengthy calculation, can be found in Appendix A.

A. Model parameters and form factors
To fully define the model it is necessary to provide the values of the unknown parameters and to specify the shape of the form factors entering the nonlocal fermion currents. There are six parameters, namely, the current quark mass m c and the dimensionful coupling constants G S , G V , G 0 , G 5 and κ. Regarding the form factors, as stated in the Introduction, we will take into account the results obtained in lattice QCD for the momentum dependence of the mass and WFR in the quark propagator. Therefore, following Ref. [26], we write the effective mass m(p) as where α m is a mass parameter defined by the normalization condition f m (0) = 1. Since LQCD calculations involve various current quark masses, we have chosen to take as input the shape of the (normalized) function f m (p 2 ), taking LQCD results in the limit of low m c and smallest lattice spacing. Considering the LQCD analysis in Ref. [26], we parameterize this function by with α = 3/2. On the other hand, for the wave function renormalization we use the parametrization [11,13] where f z (p 2 ) = 1 It is found that LQCD results favor a relatively low value for the exponent β, therefore we take here β = 5/2, which is the smallest exponent compatible with the ultraviolet convergence of the gap equations (12). As required by dimensional analysis and Lorentz invariance, the functions f m (p 2 ) and f z (p 2 ) carry dimensionful parameters Λ 0 and Λ 1 , which represent effective cutoff momenta in the corresponding channels. Thus, we will use here the above functional forms for the form factors, taking Λ 0 and Λ 1 as two further free parameters of the model. Regarding the parameters α m and α z introduced in Eqs. (61) and (63), from Eqs. (11) it is seen that they are related to the mean field values of the scalar fields by hence, for a given set of model parameters, they can be obtained by solving the gap equations (12).
The model also includes the form factors h(p), h 0 (p) and h 5 (p), introduced through the vector and axial vector current-current interactions. For definiteness and simplicity we will assume the effective behavior of quark interactions to be similar in the J = 0 and J = 1 channels, therefore we will take for h(p) the same form as g(p). Regarding the vector-isoscalar sector, as it is usually done we assume approximate degeneracy with the vector-isovector part, hence we take h(p) ≃ h 0 (p).
The axial vector-isoscalar sector can be studied separately, since it decouples from the rest of the Lagrangian. Here we will just take h 5 (p) = h(p) in order to get an estimation for the constant G 5 from phenomenology.
Given the form factor shapes, in order to study the phenomenology we have to determine the values of the model parameters (current quark mass, coupling constants and effective cutoff momenta). To do this, we first carry out a fit to lattice results for the functions f m (p 2 ) and z(p), from which we obtain the values of the cutoffs Λ 0 and Λ 1 , as well as the parameter α z . The latter will be used, together with five phenomenological quantities, as input to determine the remaining six free model parameters. From the LQCD results quoted in Ref. [4] we obtain with χ 2 /dof = 1.17 and χ 2 /dof = 0.25 for the fits to f m (p 2 ) and z(p) data, respectively. The fits have been carried out considering lattice values up to 2.5 GeV. Both the data and the fitting curves for f m (p 2 ) and z(p) are shown in Fig. 5. In the case of z(p), it is seen that the fit leads to somewhat large values of z(p) at low momenta in comparison with lattice points. We notice, however, that errors in this region are relatively large, and in addition these points are the most sensitive to changes in lattice spacing and/or sea quark masses [4]. Once the form factor shapes have been fixed, one can set the model parameters so as to reproduce the empirical values of some selected observables. As stated, we take from the fit the values of Λ 0 and Λ 1 and then we determine the values of the parameters m c , G V , G S , G 0 , G 5 and κ from six input quantities. These have been chosen to be the fitted value of α z together with the empirical values of the pion weak decay constant f π and the masses of the π, ρ, ω and f 1 mesons. From our numerical analysis we find that there is a set of parameters that allows us to properly reproduce these empirical values. The corresponding results are quoted in Table I.   Depending on the value of the three-momentum p these cuts can occasionally cross the real axis, and have to be taken into account through a proper deformation of the integration path. Details of the calculations are given in Appendix B.
From Table I we find a ratio G S /G V ∼ 1.5, which is in agreement with standard NJL model parametrizations [2]. Concerning the value of G 0 , it is necessary to take into account that we are working within a two-flavor model, therefore effects of strange quark bound states are not explicitly considered. Our determination of G 0 would be valid only in the case of "ideal mixing" between SU(3) f singlet and octet I = 0 states, which means taking the ω as an approximate SU(3) f octet state, and the φ meson as an approximately puress state. In the case of the f 1 axial vector meson there is an additional problem, which is common to various quark models. Indeed, models that do not include an explicit mechanism of confinement usually have difficulties for describing meson resonances, since there is a threshold above which the meson mass becomes large enough to allow the decay of the meson into two quarks. This threshold is typically of the order of 2m(0), therefore models that lead to constituent masses larger than about 400 MeV (as occurs in our case) can avoid this problem for low mass resonances like the ρ meson [27]. Other possible approaches are e.g. the extension of G M (−s) functions to the complex plane [28] or the search for a peak in the meson spectral function [29]. Mathematically, in our model the onset of the unphysical qq channel corresponds to the fact that in the integrals of the form of e.g. Eq. (17) there is a "pinch point" at which both functions D(q + ) and D(q − ) in the integrand are equal to zero (i.e. both constituent quarks are simultaneously on shell). For the parameters in Table I Our numerical results for various observables are summarized in Table II (we have not included here the quantities taken as phenomenological inputs, namely m π , f π , m ρ , m ω and m f 1 ). From the table it is seen that the predictions of the model for the π 0 → γγ, ρ → e + e − and ρ → ππ decay rates are in good agreement with experiments, being compatible with the empirical values [30] within an accuracy of less than 10%. We can also obtain a prediction for the width Γ(ω → e + e − ), which is found to be about 0.8 keV, somewhat larger than the experimental value 0.60 ± 0.02 keV [30].
However, as discussed above, our result might become modified after the inclusion of strangeness degrees of freedom owing to the ω − φ mixing. Regarding the σ − σ ′ sector, we obtain a physical state with a mass of about 680 MeV, which can be identified with the observed σ meson resonance (the mass of which is rather uncertain), while for the state σ ′ we find that the function G σ ′ (−s) grows monotonically with s, indicating that this state does not represent a physical meson (a more detailed discussion on the σ ′ state in this type of models can be found in Ref. [13]). In the case of the a 1 vector mesons we find that the function G  Finally, in Table III we quote our results for mean field values of scalar fields, chiral quark condensates and effective quark-meson couplings. It is seen that the model leads to a zero-momentum effective quark mass m(0) = (m c +σ 1 )/(1 −σ 2 ) ≃ 400 MeV, somewhat larger than the value of 311 MeV obtained in Ref. [13] for a nlNJL model without vector meson degrees of freedom.
For comparison, notice that standard NJL model parametrizations lead to values of constituent (momentum-independent) quark masses around 350 MeV [2]. Concerning the chiral quark condensates, our results are relatively large in comparison with usual phenomenological estimations and lattice calculations, which lead to condensates in the range of (−240 MeV) 3 to (−320 MeV) 3 [31].
In addition, when determining the model parameters we have found a relatively low value for the current quark mass, namely m c = 1.59 MeV, in comparison with lattice estimates that lead to m c ≃ 3.4 ± 0.25 MeV in the isospin limit [30]. The results for these quantities in nlNJL models are in fact strongly dependent on the form factor shapes, as it is found in Refs. [13,20,21], where two-and three-flavor nonlocal models (which do not include the vector meson sector) are considered. As discussed in those articles, one has to take into account that both m c and qq are scale-dependent quantities, and our fit has been carried out using lattice data that correspond to a renormalization scale µ = 3 GeV, somewhat larger than the usual scale of 2 GeV. To get rid of the scale dependence one can look at the product − qq m c , for which we get, within our parametrization, a result of about 8.12 × 10 −5 GeV 4 . This is in good agreement with the value arising from the Gell-Mann-Oakes-Renner relation at the leading order in the chiral expansion, namely − qq m c = f 2 π m 2 π /2 ≃ 8.21 × 10 −5 GeV 4 . Finally, for completeness we include in Table III the values obtained for the effective quark-meson couplings g πqq and g ρqq .  It is worth mentioning that analytical expressions for the pion mass and decay constants obtained in previous works have been revisited in order to take into account π−a 1 mixing.
On the phenomenological side, the fit of nonlocal form factors to lattice QCD results for effective quark propagators provides a more natural and realistic way to regularize the model in comparison with the standard NJL approach. The remaining unknown parameters, namely the current quark mass and the current-current coupling constants, can be determined from some input observables.
Here we have chosen to take as inputs the measured values of the pion decay constant and a set of meson masses. From the numerical evaluation of the analytical expressions we find that the model is able to properly reproduce the empirical values of these observables, and leads to phenomenologically acceptable values for other scalar and vector meson masses and decay widths.
To conclude, let us state that the inclusion of the axial and vector meson sector offers a more complete picture of hadron phenomenology in the framework of nonlocal quark models, and its effects can be important for the analysis of hadronic observables such as the pion electromagnetic form factor and the vector and axial vector form factors for pion radiative decays. It is also worth it to extend the study of ρ meson properties to finite-temperature systems, given its importance for the study of heavy ion collisions. In addition, for the case of hadronic systems at finite chemical potential it is expected that vector interactions lead to a nonzero condensate in the J = 1, I = 0 channel, which can be important for the study of the QCD phase diagram [32] and the physics of compact objects [33]. We expect to report on these issues in forthcoming articles.
To perform the calculations we choose, as usual, the 4th axis in the direction of the external momentum. Thus one has p µ = (iM, 0), and I(p 2 ) can be reduced to a double integral in q 4 and | q|. Since the functions F (q + , q − , p) are symmetric under the exchange q + ↔ q − , it is easy to see that F (q + , q − , p) = F (q + * , q − * , p), which ensures the reality of I(q 2 ). Now let us take | q| fixed, and consider the analytical structure of the integrand in the complex q 4 plane. It is immediately seen that we will find a pair of branch cuts in this plane arising from the function f m (q 2 ), and other pairs of cuts will appear from the occurrences of f m [(q + ) 2 ] and f m [(q − ) 2 ], respectively. In the case of f m (q 2 ) = f m (q 2 4 + | q| 2 ), the cuts are given by Re(q 4 ) = 0, |Im(q 4 )| > | q|, hence they never cross the real q 4 axis, along which the integral is to be performed. On the other hand, for The treatment of these cuts is a matter of prescription. In fact, after taking the form factors from LQCD calculations in Euclidean space, one could turn back to Minkowski space through a Wick rotation. Then one would find that the cuts are located along the integration axis, and to evaluate the integrals they have to be moved away according to some recipe. Here we will adopt the prescription of translating the arguments of f m (s) according to while f m (q 2 ) is kept unchanged. In this way, branch cuts do not overlap and the property The corresponding curves in the complex plane q 4 are sketched in Fig. 6, where we have distinguished two situations in which | q| > M/2 (Fig. 6a) and | q| < M/2 (Fig. 6b). Branch cuts corresponding to the functions f m [(q + ) 2 ], f m [(q − ) 2 ] and f m (q 2 ) have been represented with dashed, dotted and dashed-dotted lines, respectively. If | q| > M/2, as it is shown in Fig. 6a, the cuts do not cross the integration axis, thus there is no extra contribution to the loop integral. On the contrary, for | q| < M/2 two branch cuts cross from one half-plane to the other one, passing through the real q 4 axis. Since the integral over q 4 has to be ultimately equivalent to an integral over the Minkowski momentum q 0 , obtained through the corresponding Wick rotation, the integration contour along q 4 should be deformed in order to subtract the contribution of the crossing pieces, which are represented with solid lines in Fig. 6b. A similar procedure has to be followed when poles of the integrand cross the integration axis at some value of | q|; in that case the contributions resulting from the deformation of the q 4 integration contour can be obtained by calculating the residues of the poles, according to Cauchy's theorem. The need to add cut or pole contributions to the loop integrals becomes evident by looking at relatively simple integrals as those appearing in the gap equations (12): if one carries out a translation of the loop momentum p → p ′ = p + r, with r 2 = −M 2 , for fixed | p ′ | there will be branch cuts in the complex plane p ′ 4 that cross from the upper half-plane to the lower one (or vice versa). In addition, in general the integrand will have poles that for large enough values of M cross the real p ′ 4 axis at some value of | p ′ |. From Cauchy's theorem it is easy to see that the corresponding contributions have to be subtracted if one requires the loop integral to be invariant under the translation. In practice the contributions from the cuts can be obtained by carrying out integrations in the q 4 plane along adequate contours that enclose the crossing pieces, letting then ε → 0. Owing to