Thermo-magnetic spectral properties of neutral mesons in vector and axial-vector channels using NJL model

In this work the neutral meson properties have been investigated in the presence of thermo-magnetic background using two-flavor Nambu–Jona-Lasinio model. Mass, spectral function and dispersion relations are obtained in the scalar (σ) and pseudo-scalar (π0) channels as well as in the vector (ρ0) and axial vector (a0 1 ) channels. The general Lorentz structures for the vector and axial-vector meson polarization functions have been considered in detail. The ultra-violet divergences appearing in this work have been regularized using a mixed regularization technique where the gamma functions arising in dimensional regularization are replaced with incomplete gamma functions as usually done in the proper time regularization procedure. The meson spectral functions obtained in the presence of magnetic field possess non-trivial oscillatory structure. Similar to the scalar and pseudo-scalar channel, the spectral functions for each of the modes of ρ0 are observed to overlap with the corresponding modes of its chiral partner a 1 mesons in the chiral symmetry restored phase. We observe discontinuities in the masses of all the mesonic excitations for non-zero external magnetic field.


I. INTRODUCTION
Based on a considerable amount of research regarding the generation of magnetic fields in non-central heavy ion collision(HIC), there exists a growing consensus that an extremely strong transient magnetic field of the order of ∼ 10 18 G or larger can be produced at RHIC and LHC [1][2][3][4][5][6][7][8]. Being comparable to the energy scale of strong interaction, though short lived, the produced field can impart significant modifications in the properties of strongly interacting matter [9][10][11][12][13][14] resulting in a plethora of novel phenomena like chiral magnetic effect [2,15,16], magnetic catalysis [17][18][19][20], inverse magnetic catalysis [21][22][23] electromagnetically induced superconductivity and superfluidity [24,25] and so on. The tools and techniques developed for studying such magnetic modifications in HIC experiments also bear significant importance for their applicability in many different physical scenarios where strong magnetic field can be realized. For example, in the early universe during the electroweak phase transition, magnetic fields as high as ∼ 10 23 G might have been produced. Also, in case of magnetars surface magnetic field is of the order of ∼ 10 15 G. In the interior, the field intensity is even higher reaching up to ∼ 10 18 G. Such low temperature and high density extreme states are expected to be explored in CBM experiment at FAIR. On theoretical grounds, at lower temperatures, usual field theoretic approach of studying quantum chromodynamics (QCD) is not feasible due to the confining nature of strong interaction that severely restricts the applicability of perturbative analysis. In this scenario, an alternative to the non-perturbative Lattice QCD approach is provided by the QCD inspired effective models. The modification of such effective descriptions in presence of external magnetic field has gained significant research interests in recent times [26]. One such model is the Nambu-Jona-Lasinio model [27][28][29][30][31][32] which has been widely used in the studies of chiral symmetry breaking as well as meson properties in presence of thermo-magnetic background [33][34][35][36].
In the context of studying the mesonic properties in presence of magnetic field in NJL model, often the lightest mesons σ and π are considered [37][38][39][40][41][42]. In some studies diquarks are also included [43]. ρ meson properties have been discussed in Ref. [44][45][46]. In Ref. [44], it is observed that at vanishing magnetic field, there exists a temperature when ρ mass coincides with twice the mass of the constituent quark and beyond that temperature no solution for the ρ meson mass exists which is described as the ρ melting. Even at finite magnetic field the melting persists and two different melting temperatures are observed corresponding to the charged and the neutral ρ. Comparison with ρ 0 meson suggests that melting of ρ ± occurs at lower temperature in the presence of magnetic field. For example, in case of charged ρ, no solution exists beyond temperature 169 MeV for eB around 0.2 GeV 2 . However, similar analysis on ρ ± in Ritus formalism [46] does find non vanishing mass for charged ρ even at much higher values of temperature for similar strength of the background magnetic field ( see for example Fig.4 of Ref. [46]). The apparent ambiguity thus demands investigation of the properties of neutral ρ meson in thermo-magnetic background which essentially will be an extension of the study presented in Ref. [46].
On a different note, one of the significant features of studying the meson properties is that at temperatures higher than a critical value, the masses of the chiral partners become degenerate. This degeneracy in the meson mass spectrum serves as an important signature of the chiral symmetry restoration. Therefore, the restoration of chiral symmetry in the vector channel can be shown explicitly only when one includes the a 1 meson along with ρ which is missing in the studies of ρ mesons discussed earlier [44,46]. It should be mentioned here that in order to investigate the vector and pseudo-vector channel, proper incorporation of the general structure of meson self-energy is required. The general Lorentz structure for the ρ meson in presence of thermo-magnetic medium has been recently reported in Ref. [47]. One may note that the Lorentz structure of ρ meson polarization function has not been considered in Ref. [44].
In this work the neutral meson properties in scalar (σ) and pseudo-scalar (π 0 ) channels as well as vector (ρ 0 ) and axial vector (a 0 1 ) channels have been investigated in the framework of two-flavor NJL model in presence of constant background magnetic field. The detailed general structure for the vector and axial-vector meson polarization tensor have been considered. The Schwinger scheme has been implemented in the evaluation of the polarization tensors. However, as only the neutral mesons are considered, the Schwinger phase vanishes and Schwinger and Ritus formalisms are expected to provide identical results [41]. It should be mentioned here that being an effective description of QCD at low energy regime, NJL model is non-renormalizable and requires a regularization prescription. The most commonly used regularization technique is to use a three momentum cut-off which acts as a parameter of the theory and can be fixed to reproduce some well-known phenomenological quantities, for example the pion-decay constant and the condensate value. However, to obtain the general structure of the self energy in a consistent way, we take recourse to dimensional regularization(DR) technique. Now, the ultra-violet divergences in dimensional regularization prescription occurs as pole of gamma function. In that procedure, one extra parameter arises which is to be simultaneously fitted to reproduce the phenomenological quantities. Detailed description regarding the fitting procedure can be found in Refs. [48][49][50]. However, in this work, to obtain the finite contribution, the gamma functions arising from DR are replaced with incomplete gamma functions. We refer this replacement procedure as incomplete gamma regularization. As a reward, though the number of parameter set remains identical to that of usual regularization procedures, in this scheme, the general Lorentz structure for vector and axial vector polarization functions can be obtained systematically. The regularization scheme has been used to obtain the neutral meson properties like mass, spectral function and dispersion relations. Non-trivial mass jump is observed in the spectrum for each of the modes in the vector and axial vector channel which bears similarity with earlier studies of pions in presence of magnetic field [41,42].
The article is organized as follows. Sec. II describes the constituent quark mass and the dressed quark propagators in the real time formalism of thermal field theory whereas the gap equations and general structure are described in Sec. III. In both the sections, vacuum, thermal and thermo-magnetic cases are considered in separate subsections. The main results for the real and imaginary parts of the meson polarization functions are listed in Sec. IV. Sec. V describes the regularization procedure used in this work. All the numerical results are presented in Sec. VI followed by a brief summary in Sec. VII. Some of the relevant calculational details are provided in the appendices.

II. THE CONSTITUENT QUARK MASS AND THE DRESSED QUARK PROPAGATOR
The standard expression of the two-flavour NJL Lagrangian is where, ψ = u d is the quark isospin flavour doublet with u and d being the up and down quark fields respectively. Each of the up and down quark fields are [4 × 1] ⊗ [3 × 1] matrices corresponding to their orientation in Dirac and color spaces. In Eq. (1), g s and g v are respectively the coupling constants in the spin-0 and spin-1 channels for the four point contact interactions among the quark fields and m is the current quark mass which is assumed to be equal for the up and down quarks to ensure isospin symmetry. In the NJL model the constituent quark mass is dynamically generated as a consequence of the spontaneous breaking of chiral symmetry. In the following subsections, we briefly introduce the formalism required to obtain the constituent quark mass and the dressed quark propagator for three different cases separately: We first consider the pure vacuum case for which the temperature is zero and external magnetic field is switched off. The dressed quark propagator S ′ (q) is calculated from the Dyson-Schwinger equation where, S(q, m) = −( ✁ q+m) q 2 −m 2 +iǫ ⊗ 1 Flavor ⊗ 1 Color is the free quark propagator and Σ is the one-loop self energy of quark. In the Mean Field Approximation(MFA), the quark self energy becomes diagonal in Dirac, color and flavor spaces as This enables one to solve Eq. (2) trivially to get the complete propagator as where is the 'constituent quark mass'. The above equation is the well known 'gap equation'. Our next task is to calculate the quantity Σ MFA . Applying Feynman rules to Fig. 1, we get the one-loop self energy of quark in the MFA as It is to be noted that, the loop particle in the self energy is dressed. In the above equation, the subscripts c, d and f in the Tr correspond to the traces taken over color, flavor and Dirac spaces respectively. Also note that, the quark self energy is a function of M itself (since the loop particle is dressed) so that Eq. (5) has to be solved self-consistently to calculate M.
Let us now explicitly evaluate the quantity Σ Pure-Vac MFA . Substituting Eq. (4) into Eq. (6), we get where, N c = 3 and N f = 2 are the number of colors and flavors respectively. The momentum integral in the above equation is ultra-violet (UV) divergent. The NJL model, being a non-renormalizable theory, requires a proper regularization scheme. There exists many such UV-regulators in the literature such as three-momentum cutoff, Euclidean four-momentum cutoff, Pauli-Villars, proper time and so on. The mostly used regulator is the momentum cutoff which breaks the Lorentz invariance and usually every symmetry of the theory. It will be demonstrated later in Sec. V that, the momentum cutoff regulator (or any other regulator which breaks Lorentz invariance) is not useful to study the vector meson ρ in the NJL model. In this work, we will use 'dimensional regularization' as our UV-regulator which respects all the symmetries of the theory. Going to d-dimension, Eq. (7) becomes where λ is a scale of dimension GeV 2 which has been introduced to keep the overall dimension of the equation consistent.
Performing the momentum integral in the above equation, we get where ε = 2 − d/2. It is to be noted that, the UV-divergence has been isolated as the pole of the Gamma function since Γ(z) has simple poles at z = 0, −1, −2, ..... The regularization procedure of the above divergent quantity will be discussed in Sec. V. The above equation has the following expansion about ε = 0 which will be used later.
We now turn on the temperature and consider the case T 0, B = 0. To include the effect of finite temperature, we will use the Real Time Formalism (RTF) of finite temperature field theory [51,52]. In the RTF, all the two point correlation functions including self energies and propagators become 2 × 2 matrices (will be denoted by boldface letters) in thermal space. As a result the Dyson-Schwinger equation generalizes to a matrix equation in thermal space where each of the quantities is a 2 × 2 matrix. In the above equation, S is the free thermal quark propagator given by In the above equation, the diagonalizing matrix V is given by, where, u µ is the four velocity of the thermal medium. In the Local Rest Frame (LRF), one has u µ LRF ≡ (1, ì 0). In the above equations, Θ(x) is the unit step function and f (x) = e x/T + 1 −1 is the Fermi-Dirac thermal distribution function for the quarks.
It is well known that, the complete thermal propagator matrix S ′ and the thermal self energy matrix Σ are diagonalized by V and V −1 respectively. Thus Eq. (11) boils down to an algebraic equation in thermal space as where S ′ and Σ are respectively the 11-component of the matrices VS ′ V and V −1 ΣV −1 . As before, in the MFA, the Σ is diagonal in Dirac, color and flavor spaces Σ = Σ MFA 1 Dirac ⊗ 1 Flavor ⊗ 1 Color so that Eq. (16) can be trivially solved to obtain S ′ (q, m) = S(q, M(T )) where the thermal constituent mass is given by It is easy to check that, ReΣ = ReΣ 11 , so that the knowledge of S ′ 11 is sufficient to calculate the quantity ReΣ. The explicit form of S ′ 11 is given by where Let us now evaluate the thermal self energy function Σ MFA whose real part is obtained by replacing the vacuum complete propagator on the RHS of Eq. (6) by S ′ → S ′ 11 as Substituting Eq. (19) into the above equation, we get after some simplification where ω k = ì k 2 + M 2 .
C. CASE-III: T 0, B 0 Finally, we consider the case of finite temperature and non-zero external magnetic field. In this case, the complete thermomagnetic quark propagator S ′ B satisfies the generalized Dyson-Schwinger equation where Σ B is the thermo-magnetic quark one-loop self energy matrix and S B is the free thermo-magnetic quark propagator. Analogous to Eq. (12), S B can be written explicitly as where, in which, S u and S d are respectively the Schwinger proper-time propagator for up and down quarks. They can be expressed as a sum over discrete Landau levels as where, α and with e f being the electric charge of flavour f i.e. e u = 2 3 e and e u = − 1 3 e; e is the charge of a proton. In the above equation, sign (x) = Θ(x) − Θ(−x); L a l (z) is the generalized Laguerre polynomial with the convention L a −1 = 0. The external magnetic field being in the positive z-direction, the metric tensor can be decomposed as g µν = g µν + g µν ⊥ where g µν = diag(1, 0, 0, −1) and g µν ⊥ = diag(0, −1, −1, 0) so that the parallel and perpendicular four vectors are defined as q µ = g µν q ν and q µ ⊥ = g µν ⊥ q ν . Similar to the thermal case, the Dyson-Schwinger equation in thermo-magnetic medium can also represented in diagonal form as Following the MFA, the Σ is diagonal in Dirac, color and flavour spaces, Σ = Σ MFA 1 Dirac ⊗ 1 Flavour ⊗ 1 Colour so that Eq. (28) can be trivially solved to obtain S ′ (q, m) = S B (q, M(T, B)) where the thermo-magnetic constituent quark mass is given by As before, because of the fact ReΣ = ReΣ 11 B , the knowledge of S ′ B11 is sufficient to calculate the quantity ReΣ. The explicit form of S ′ B11 is given by where Let us now evaluate the thero-magnetic self energy function Σ MFA whose real part is obtained by replacing the 11-component of the complete thermal propagator on the RHS of Eq. (20) Substituting Eq. (30) into the above equation, we get after some simplification (see Appendix A for details) where, ReΣ B-Vac MFA (M, B) is the real part of the magnetic field dependent vacuum contribution to the quark self energy and can be read off from Eq. (A13) as The temperature as well as magnetic field dependent contribution to the self energy, ReΣ B-Med MFA (M, B, T ), can be obtained from Eq. (A7) as It is interesting to note that in Eq. (34), the divergent pure vacuum contribution ReΣ MFA (M) has completely been decoupled from the magnetic field and temperature dependent parts. One can notice from Eq. (35) that, the quantity Σ B-Vac MFA (M, B) is finite and thus the external magnetic field does not produce any additional divergences.
The formalism described in Appendix A to untangle the divergent pure vacuum contribution of the one-loop self energy graph is closely related to the Magnetic Field Independent Regularization (MFIR) scheme as developed in Ref. [50]. However, the methodology we have adopted is slightly different from that of MFIR. In our case, we have performed a dimensional regularization to the d d k integral which leads to Hurwitz zeta function as a function of the dimension d. An expansion of the Hurwitz zeta about its pole leads to the disentanglement of the pure vacuum part. On the other hand, in MFIR scheme, one does not change the space time dimension; rather one adds and subtracts the pure vacuum part. Then, using an integral representation of the Hurwitz zeta function, the vacuum subtracted self energy is written as an integral over Schwinger proper-time parameter . Finally the proper-time integral is performed to get the vacuum subtracted finite magnetic field dependent self-energy. Regardless of the methodology used, the two procedures lead to similar results. Specifically, the expressions in Eqs. (33)- (36) are identical to the ones obtained in Refs. [37,50].

III. MESON PROPAGATORS IN RANDOM PHASE APPROXIMATION IN THE NJL MODEL
Since mesons are the bound state of quarks and anti-quarks, their propagation can be studied from the scattering of quarks in different channels using the Bethe-Salpeter approach [29]. On the other hand, as discussed in Ref. [45,46], the meson propagators can also be recast into the form of Dyson-Schwinger equations in the random phase approximation (RPA). Let us first consider the situation at vacuum (i.e. T = 0 and B = 0). In the scalar and pseudo-scalar channel, the π and σ meson propagators D ′ h satisfy the following Dyson-Schwinger equation where D = (−2g s ) are the bare propagators and Π h are the one-loop polarization functions. The corresponding expression of the meson propagators D ′µν h in the vector (ρ) and pseudo-vector (a 1 ) channels are given by where D µν = (2g v g µν ) are the bare propagators and Π µν H are the one-loop polarization functions for the ρ and a 1 mesons. As already discussed in Sec. II, at finite temperature, all the real time two point correlation functions become 2 × 2 matrices in thermal space and will be denoted by boldface letters. Thus, at finite temperature, Eqs. (37) and (38) generalize to However, each term of the above equations can be diagonalized to express them in terms of analytic functions [52] (will be denoted by bars) which in turn diagonalizes the Dyson-Schwinger equation making it an algebric equation in thermal space as In presence of both the finite temperature and external magnetic field, the generalization of Eqs. (39) and (40) is so that, the the corresponding thermo-magnetic analytic functions denoted by a double-bar satisfy, Our next task is to solve the Dyson-Schwinger equations in order to express the complete meson propagators in terms of the polarization functions. It is trivial to solve Eqs. (37), (41) and (45) for the π and σ mesons as However, for the ρ and a 1 channels, additional complications arise because of the Lorentz indices in Eqs. (38), (42) and (46). It is useful to decompose the polarization function and the complete propagator in terms of orthogonal tensor basis (constructed using the available vectors and tensors). This will enable one to solve the corresponding Dyson-Schwinger equation in a covariant way. We will discuss this in the following subsections.

A. GENERAL LORENTZ STRUCTURE OF THE SPIN-1 POLARIZATION FUNCTION
In order to decompose the polarization function into a suitable Lorentz basis, we use the fact that the polarization function is symmetric in its two Lorentz indices. We start with the simplest case of vacuum (i.e. T = 0, B = 0). The available quantities to construct a tensor basis are the momentum of the meson q µ and the metric tensor g µν . Only two basis tensors can be constructed which are the following It can easily be checked that, P µν i with i = 1, 2 satisfies all the properties of projection tensors i.e.
The vacuum polarization function Π µν H in this basis can be written as where the form factors Π H i are obtained using Eq. (51) as Note that, the form factors can be expressed in terms of the Lorentz invariants that can be formed by contracting Π µν H with the available tensors and vectors. In this case, we have g µν and q µ so that the form factors can be expressed in terms of g µν Π µν H and q µ q ν Π µν H . See Appendix C for details. Let us now consider the case of finite temperature only (i.e. T 0 and B = 0). Apart from q µ and g µν , in this case we have an additional four-vector u µ . Thus one can choose the following four tensors as the basis is a vector orthogonal to q µ . Similar to the vacuum case, one can verify that, the above tensors qualify to be the orthogonal projection tensors as they satisfy where the angular bracket is the shorthand notation for A µ B ν = A µ B ν / √ A 2 B 2 . Now, the analytic thermal polarization function Π µν H can be expanded in the above basis as where the form factors Π H i are obtained using Eq. (60) as Note that, the form factors can be expressed in terms of the Lorentz invariants that can be formed by contracting Π µν H with the available tensors and vectors. In this case, we have g µν , q µ and u µ so that the form factors can be expressed in terms of Significant care has to be taken while considering the special case of ì q = ì 0 [47,52]. To see this, let us write ì q = | ì q|n wheren is the unit vector in the direction of ì q. In the limit of | ì q| → 0, we have implying that the above components of the polarization tensors depends of the direction of ì q even if | ì q| = 0. This ambiguity is rectified by imposing additional constraints on the form factors Π i as Finally, we consider the general case of both finite temperature as well as finite external magnetic field. In this case, another four vector b µ = 1 B G µν u ν appears which specify the direction of the external magnetic field in the LRF where G µν = 1 2 ǫ µναβ F αβ is the dual of the field tensor F αβ (we have used ǫ 0123 = 1). In the LRF, we have b µ LRF ≡ (0, 0, 0, 1). Thus using q µ , u µ , b µ and g µν , we can construct the following seven orthogonal tensors: It can be shown that, the tensors P g µν g αβ Now, the analytic thermo-magnetic polarization function Π µν H can be expanded in the above basis as where the form factors Π H i are obtained using Eq. (75) as As before, the form factors can be expressed in terms of the Lorentz invariants that can be formed by contracting Π µν H with the available tensors and vectors. In this case, we have g µν , q µ , u µ and b µ so that the form factors can be expressed in terms of seven Similar to the thermal case, significant care has to be taken while considering the special case of ì q ⊥ = ì 0. To see this, let us write ì q ⊥ = | ì q ⊥ |n wheren is the unit vector in the direction of ì q ⊥ . In the limit of | ì q ⊥ | → 0, we have which implies that the above components of the thermo-magnetic polarization tensors depends of the direction of ì q ⊥ even if | ì q ⊥ | = 0. This ambiguity is rectified by imposing additional constraints on the form factors Π i as

B. SOLUTION OF THE DYSON-SCHWINGER EQUATION AND COMPLETE SPIN-1 PROPAGATORS
Having obtained the general Lorentz structure of the polarization functions in the previous subsection, we can now solve the Dyson-Schwinger Eqs. (38), (42) and (46) in order to calculate the complete propagators for ρ and a 1 mesons.
Let us start with solving Eq. (38). We first write where the form factors D H i are to be determined. Rewriting Eq. (38) as and making use of (50), one obtains the form factors of the complete propagator as Let us now proceed to obtain the complete thermal propagator by solving Eq. (42). Expressing the complete propagator in the orthogonal tensor basis as where the form factors D H i are to be determined. Rewriting Eq. (42) as and making use of i along with Eq. (59), one obtains the form factors of the complete thermal propagator as Finally we calculate complete thermo-magnetic propagator by solving Eq. (46). Expanding the complete propagator in the orthogonal tensor basis as where the form factors D H i are to be determined. Rewriting Eq. (46) as and making use of i along with Eq. (74), one obtains the form factors of the complete thermo-magnetic propagator as where

IV. POLARIZATION FUNCTIONS OF THE MESONS
In this section, we will explicitly calculate the polarization functions in various channels. In the current work, we only include the charge-neutral mesons i.e. π 0 , σ, ρ 0 and a 0 1 . Thus by π, ρ and a 1 we will mean π 0 , ρ 0 and a 0 1 . We start with the well known expression for the vacuum polarization functions (at T = 0 and B = 0) of the charge-neutral mesons where τ 3 = 1 0 0 −1 is the third Pauli isospin matrix and S ′ (q, m) is defined in Eq. (4). Similar to the case of quark self energy calculation, we will use dimensional regularization for the evaluation of the above pure-vacuum polarization functions. The calculation has been briefly sketched in Appendix. D and the final result can be read off from Eqs. (D8)-(D11) as As finite temperature, the analytic thermal polarization functions Π h (q) and Π H (q) are related to the 11-components of respective thermal polarization matrices Π 11 h (q) and Π 11 H (q) via relations [51,52] ReΠ h (q) = ReΠ 11 ImΠ h (q) = sign(q 0 ) tanh q 0 2T ImΠ 11 h (q) and ImΠ µν H (q) = sign(q 0 ) tanh Now, the 11-components of the thermal polarization functions are obtained by replacing the vacuum propagators on the RHS of Eqs. (101)-(104) by S ′ → S ′ 11 where S ′ 11 is defined in Eq. (19). Therefore, Substituting S ′ 11 from Eq. (19) into the above equation and making use of Eqs. (109) and (110), we get after some simplifications the real parts of the polarization functions as and imaginary parts as where N h (q, k) and N Substituting S ′′ 11 from Eq. (30) into the above equation and making use of analogous relations to Eqs. (109) and (110), we will obtain the real and imaginary parts of the analytic thermo-magnetic polarization functions. For the simplicity in analytic calculations, we take q ⊥ = 0 for which the corresponding calculations are provided in Appendix E and below we only give the final expressions. From Eqs. (E18), (E19) and (E31)-(E36), we get the real parts of the analytic thermo-magnetic polarization functions as where the novel magnetic field dependent vacuum contributions are The imaginary parts are to be read off from Eqs. (E15) and (E16) as where N lnf h (q, k) and N lnf µν H (q, k) are defined in Eqs. (E11) and (E12). It may be emphasized that though the present work uses Real Time version of thermal field theory, use of the more popular imaginary time formalism leads to the same expressions. For example, the expression of the thermo-magnetic quark self energies or the polarization functions of π 0 and σ obtained here are identical to the ones obtained in Refs [29,37,44] earlier using the ITF.

V. REGULARIZATION PROCEDURE FOR THE NJL MODEL
As already mentioned in the previous sections, that the NJL model requires a proper regularization procedure. Using the dimensional regularization technique, we have been able to isolate the UV-divergences as the pole of Gamma functions in Eqs. (10) and (105)-(108). Now, in order to obtain a finite contributions from these equations, we first note that the integral representation of the Gamma functions can be written as where, γ(z, r) is the lower incomplete gamma function and Γ(z, r) is the (upper) incomplete Gamma function. In the evaluation of loop diagrams in the NJL model using Schwinger proper-time method, one often encounters integrals which can be written in terms of Γ functions with negative integer argument. Clearly, those are divergent quantities and need to be regulated. One possible way is to introduce proper-time regulator where the lower incomplete gamma function containing the divergence is discarded and only the Γ(z, r) part is retained (see for example Eq. (3.15) in [29] which is the proper-time regularized version of Eq. (3.13) there in). Following the similar procedure, in our regularization scheme, the divergent Gamma functions obtained from dimensional regularization are replaced with the incomplete Gamma function i.e.
where Λ is a scale parameter to be determined. Thus our regularization scheme is a mixed procedure where, though the dimensional regularization is used at first to obtain the consistent Lorentz structure, the divergences appeared are regulated following the proper-time regularization. After these replacements, Eqs. (10) and (105)-(108) can be simplified to and It can be notice in Eqs. (134)-(137), that if the chiral symmetry is completely restored (i.e. M = 0), then the polarization functions of σ and a 1 become identical to that of π and ρ respectively. Moreover, observing the Lorentz structure in Eq. (136), it immediately follows that, the polarization function of ρ is transverse i.e.
The reason behind this transversality is the conservation of the vector current J µ (x) = Ψ(x)γ µ Ψ which is the Noether's current corresponding to the U(1) symmetry of the NJL Lagrangian in Eq. (1). Similar arguments also hold for the Lorentz structure of the polarization function of a 1 in which the non-transverse piece is proportional to the constituent quark mass M. This is because of the non-conservation of the axial-vector current J 5µ = Ψ(x)γ µ γ 5 Ψ whose four-divergence is In the chiral limit (M = 0), the axial-vector current is conserved leading to a transverse polarization function of a 1 . It is worth mentioning that, the consistent Lorentz structure of the polarization functions of ρ and a 1 could be obtained only because we have used dimensional regularization technique which respects the Lorentz symmetry. Any other regulator such as three-momentum cutoff, Euclidean four-momentum cutoff and Schwinger proper-time regulator will spoil the Lorentz structures and Π µν ρ (q) will no longer be transverse. We now fix the parameters for the NJL model. For this we need the expression of pion decay constant ( f 2 π ) which comes out to be using the dimensional regularization. By simultaneously fitting the vacuum quark condensate and pion decay constant values as  136) and (137) will not appear automatically as appears in dimensional regularization. Moreover, the transversality condition q µ Π µν ρ (q) = 0 is not satisfied if the proper-time regularization is used.

VI. NUMERICAL RESULTS
We start this section by showing the variation of the constituent quark mass as a function of temperature for different values of external magnetic field in Fig. 2(a). As can be seen in the figure, M remains almost constant in the low temperature region. However, with further increase in temperature, the constituent quark mass decreases substantially signifying a phase transition. Throughout the whole temperature range M remains single-valued depicting the smooth crossover nature of the phase transition.
Since we are working with finite current quark mass m 0, the chiral symmetry is only partially restored. To obtain the transition temperature, one can use various susceptibilities which will be discussed in the next paragraph. For a particular value of temperature, the constituent quark mass increases with the external magnetic field as shown in Fig. 2(b).
The transition temperature corresponding to the partial restoration of chiral symmetry can be obtained from various susceptibilities. The calculation of the susceptibility −∂ M/∂T and chiral susceptibility χ = 1 external magnetic field. The position of the peak of −∂ M/∂T or χ represents the transition temperature. As can be noticed from the plots, with the increase in external magnetic field the peak of the susceptibilities moves towards higher values of temperature. Thus, in this framework, the transition temperature increases with B. This may be identified as magnetic catalysis (MC) in the NJL model where the external magnetic field catalyzes the spontaneous breaking of chiral symmetry [17][18][19][20]. Moreover, as the susceptibilities remain continuous and finite with the change in temperature, the nature of the phase transition can be inferred as smooth crossover.  We now turn our attention to the mesonic properties. We define the spectral functions of mesons as the imaginary part of the respective complete propagators. From Eq. (47), the spectral function for the π and σ mesons can be written as In Figs. 4(a)-(d), the spectral functions of π have been shown as a function of its invariant mass q 2 for different values of temperature and external magnetic field in the rest frame of π (i.e. ì q = ì 0). Let us first consider the B = 0 cases which are shown as solid-red curves in Figs. 4(a)-(d). At zero temperature, S π is a Dirac delta function at its pole mass ( q 2 = 135 MeV) along with a two-quark continuum starting at q 2 > 2M. It can be noticed from Fig. 4(b), that at T = 150 MeV, the Dirac delta function moves towards the higher invariant mass and the two-quark continuum threshold has significantly decreased which is due to the decrease in M with temperature. Yet, the delta function is well separated from the continuum revealing the fact that π is still a bound state. With further increase in temperature, as shown in Figs. 4(c)-(d), the Dirac delta function disappears and the shape of spectral function becomes a Breit-Wigner. These imply that, the pion has now become a resonant state with finite decay width. Let us now discuss the effect of external magnetic field on S π . For the lower temperature (T = 0 and 150 MeV), the Dirac delta functions move towards higher values of the invariant mass with the increase in external magnetic field. For higher values of temperature (T = 200 and 250 MeV), the spectral functions at non-zero B are observed to oscillate about the B = 0 curve and the peak of the Breit-Wigner shifts significantly towards higher invariant mass. The oscillation frequency (amplitude) is observed to be large (small) at lower values of B as compared to its higher values.
The situation is quite different in case of σ meson. In Figs. 4(e)-(h), the spectral functions of σ have been shown as a function of its invariant mass for different values of temperature and external magnetic field for ì q = ì 0. In this case, the spectral function is always Breit-Wigner shaped implying that the σ remains always a resonant excitation. As shown in Figs. 4(e)-(g), with the increase in temperature (up to T = 200 MeV), the peak of S σ moves towards lower invariant mass. However in Fig. 4(h), (at T = 250 MeV), the peak again start moving towards higher q 2 values. The effect of external magnetic field on S σ is similar to that of π showing oscillations in S σ at non-zero B about the B = 0 curve. The oscillation frequency (amplitude) follows the similar trend as described for pion.
Let us now consider the propagation of ρ and a 1 meson. Since we will be considering the special case q ⊥ = 0, we have significant simplifications of the complete propagators of ρ and a 1 . As given in Eq. (82), we have for q ⊥ = 0, Moreover, we find in our numerical calculations that Π H 5 (q 0 , ì q ⊥ = ì 0) = 0. Thus, the form factors for the complete thermomagnetic propagators in Eqs. (92)-(98) simplify to Therefore, the complete thermo-magnetic propagator from Eq. (90) becomes The second term on the RHS of the above equation containing the non-transverse tensor P µν 2 corresponds to a non-propagating mode as the corresponding form factor does not have any pole. Thus, we find three modes of propagation of ρ and a 1 mesons in the thermo-magnetic medium; two of them are found to be degenerate (corresponding to P µν 1 and P µν 4 ). This degeneracy is solely due to our special choice of q ⊥ = 0. Thus, we are left with two distinct modes for the ρ and a 1 propagations. We call them as Mode-(A) and Mode-(B) respectively. The spectral functions for these two modes are therefore defined as, In Figs. 5(a)-(p), we have presented the spectral functions of ρ and a 1 mesons as a function of their invariant mass at ì q = ì 0 for different temperature and external magnetic field. Similar to the case of σ, the ρ and a 1 are always in resonant state so that the shape of their spectral functions remains Breit-Wigner. Since we have taken ì q = ì 0 in these plots, the two modes are degenerate for B = 0 (the solid-red curves). The external magnetic field breaks this degeneracy and we find two distinct modes of ρ and a 1 propagations even in their rest frames for non-zero values of B. With the increase in temperature, the peaks of the spectral functions move towards lower values of invariant mass. Moreover, the spectral functions at non-zero external magnetic field show highly oscillatory behaviour about the B = 0 curves. Similar to the case of π and σ, we observe higher(lower) oscillation frequency (amplitude) at lower values of B.
Till now, we have taken ì q = ì 0. To see the effect of longitudinal momentum on the spectral function, we have plotted the spectral functions of the mesons as a function of invariant mass for T = 250 MeV and q ⊥ = 0 with different values of q z and external magnetic field in Figs. 6(a)-(f). First of all, it can be observed that, the spectral functions of σ and a 1 become identical to that of π and ρ respectively in all the cases as a consequence of the chiral symmetry restoration . In all the cases, the effect of increase in the q z decreases the height of spectral functions with a marginal change of their peak positions. Moreover, comparing the green-dot and violet-dash-dot curves in Figs. 6(c) and (d), it can be noticed that, a non-zero value of q z lifts the degeneracy of the two modes of ρ and a 1 at B = 0.
We now turn our attention to the study of the effect of temperature and external magnetic field on the meson masses and dispersion relations. We define the dispersion relations of the mesons as the value of ω( ì q) at which the spectral function S(q 0 = ω, ì q) has a peak (global maxima) or in other words the locus (q 0 = ω, ì q) of the peak of the spectral function gives the dispersion relations. Thus, the (effective) masses of the mesons are obtained by putting ì q = ì 0 in the dispersion relation i.e. m h, H = ω h, H ( ì q = ì 0). In Fig. 7(a), the masses of the mesons are plotted as a function of temperature at vanishing external magnetic field. Twice of the constituent quark mass is also shown for comparison. In the lower temperature region, the meson masses remains almost constant. However, m π starts increasing monotonically with temperature beyond T ≃ 150 MeV and eventually it becomes larger than 2M. On the other hand, m σ first decreases to attain a minimum after which it increases. In the whole temperature range, m σ remains always greater than 2M maintaining its resonant signature. At high temperature, the mass of π and σ merge with each other as a consequence of the chiral symmetry restoration. Similar behaviour can also be noticed for m ρ and m a 1 where both decrease with temperature followed by a merging of their masses in the chiral symmetry restored phase.
It is to be noted that, the mass/dispersion relation of the meson (or of any unstable resonance particle) can have different definition. The mass/dispersion relation can either be obtained from the locus (q 0 , ì q) of the pole of the propagator or of the peak of the spectral function. In the current work, we have used the peak of the spectral function for the definition of mass/dispersion relation. However, to check how these two differ from each other, we have plotted the masses of σ and π 0 as a function of temperature at B = 0 in Fig. 8. As can be seen from the Fig. 8, the two different definitions of mass lead to no noticeable difference. Moreover, the ratio of the masses calculated from the pole to that from the peak is exactly unity when the particle has zero decay width (for example the π 0 mass at the low temperature). Now, keeping the temperature fixed at T = 0, the variation of meson masses as a function of external magnetic field are plotted in Fig. 7(b). Frequent mass jumps are observed for the distinct modes of ρ and a 1 . In between the two successive discontinuities, the effective mass increases with eB. It can be noticed that the frequency of oscillation decreases with the external field. In other words, separation between the two successive discontinuities increases with eB. Also in case of σ mesons, the effective mass shows increasing trend between the successive discontinuities. However, only one mass jump can be seen within the plotted range of the magnetic field. Pion mass on the other hand remain continuous and is observed to decrease slowly with the external field which is consistent with Refs. [38,41]. If Fig. 9, we have shown the variation of m π as a function of temperature at two different values of external magnetic field. At lower values of temperature, mass of pions are almost independent of T . At some particular temperature, m π suffers a sudden jump (discontinuity) corresponding to Mott transition [41,42,[53][54][55]. The jump structure is in qualitative agreement with most of the studies. However, there exists differences in quantitative nature of the jump structure. For example, the amount of discontinuity obtained here is smaller in comparison to [42] which itself is different from [53] as well as [41]. One should observe that different parameter sets have been chosen in all these cases along with different regularization procedures.
Temperature dependence of m σ is shown in Fig. 10(a) at different values of the external magnetic field. At lower values of temperature, nature of m σ is dominated by its eB dependence. Because of the mass jump present at T = 0, m σ shows nonmonotonic behaviour with respect to eB variation. For example, effective mass at eB = 0.10GeV 2 is smaller than the effective mass at eB = 0.05 GeV 2 whereas the corresponding value of m σ at eB = 0.15GeV 2 remains well above the former two cases. As a result, with the increase of temperatures, when m σ decreases, crossing between fixed eB curves develops. With further increase of temperature, effective mass shows discontinuous jump structure for eB = 0.10 and 0.15 GeV 2 . This mass jump signifies the fact that even in case of sigma meson, there exist certain set of T and eB values for which no solution exists for the pole of σ propagator. The pole reappears at a higher value giving rise to a discontinuous jump. In general, this behaviour can be attributed to the oscillatory nature of the polarization function. One important feature to be noted is that at eB = 0.05 GeV 2 , the effective mass of σ does not possess any discontinuous jump within the plotted temperature range. We have also checked in our numerical calculations that, at finite temperature as well as at non-zero magnetic field, the relation m 2 σ (T, B) = m 2 π (T, B) + 4M 2 (T, B) is in agreement with Refs. [29,37,44].
In Fig. 10(b), m ρ is plotted as a function of temperature for different values of external magnetic field. The B = 0 curve is degenerate for the two modes. The degeneracy is lifted once the external magnetic field is turned on. For a given value of eB, m ρ shows a decreasing trend with temperature except at particular values where discontinuous jump occurs. The nature of the discontinuities is similar to that of m π and m σ i.e at the point of discontinuity, the solution for the pole position always jumps to higher values. Also in this case, one can observe that there exists certain magnetic fields for which no discontinuity appears within the plotted temperature range (see for example, Mode-(B) at eB = 0.15 GeV 2 ). On the other hand, for a particular temperature, m ρ is found to be oscillatory with the change in eB. In other words, the effective mass can go to higher as well lower values depending upon the external magnetic field. This is again expected from the highly oscillatory nature of the effective mass at T = 0 (shown in Fig. 7(b) ). Analogous feature is observed for the case of a 1 meson as shown in Fig. 10(d). However, in this case, the effective mass of a 0 1 can jump to lower values as well (see for example, Mode-(A) at 0.10 GeV 2 ). Finally, we concentrate This degeneracy is lifted when we take either T or B 0. Moreover, for B = 0 and T 0, the dispersion curves are identical around q z ≃ 0. The nature of the dispersion curves at different values of eB are similar to that of π and σ as they are mostly separated near q z ≃ 0 and tend to merge at high q z . The corresponding plots for the a 1 meson is shown in Figs. 12(i)-(p) and the nature of the curves are similar to that of ρ meson.

VII. SUMMARY AND CONCLUSION
In this work, the neutral meson properties such as mass, spectral function and dispersion relations have been studied in the presence of a constant background magnetic field using two-flavor Nambu-Jona-Lasinio model. The novelty of the study lies in the detailed consideration of the general Lorentz structure for the vector and axial-vector meson polarization functions, which, to the best of our knowledge, has been ignored in similar studies of vector mesons. Apart from the consideration of the modified Lorentz structure in presence of magnetic field, the Schwinger propagator expressed as sum over Landau levels has been used in the calculation of the quark self energy and meson polarization functions. For simplicity in the analytic calculation, only longitudinal mesons (q ⊥ = 0 ) are considered. To obtain the Lorentz structure of the vector and axial-vector meson systematically, we have adopted a hybrid regularization scheme where as a first step, the dimensional regularization is used to isolate the ultra-violet divergences as the poles of gamma functions. Subsequently, those gamma functions are replaced by incomplete gamma functions as usually done in the proper time regularization scheme. We call this hybrid regularization procedure as the incomplete gamma regularization (IGR). As a reward, the number of parameters remain identical to that of usual cut-off regularization procedures. We have obtained two distinct modes for ρ 0 and a 1 meson. At eB = 0 the effective mass of the modes remain degenerate, however, external magnetic field lifts the degeneracy. At temperatures above the critical temperature for chiral symmetry restoration, the spectral functions for each of the modes of ρ 0 are observed to overlap with the corresponding modes of its chiral partner a 0 1 meson for both zero and non-zero values of external magnetic field. The discontinuity in the pion mass near the Mott transition temperature is observed which is consistent with recent works [41,42]. However, in our case, the discontinuous mass jump is also observed in the effective mass of sigma meson which seems to be absent in Ref. [38](see Fig.1). Also in [41], it is mentioned that no mass jump for σ can exist in NJL model as m σ always lies above 2M. In our work too, we observe that the condition m σ > 2M is always satisfied. Thus, we conclude that this condition may not be the correct explanation of the absence of mass jump in case of m σ in [38]. In our work, discontinuous mass jumps have also been observed in different modes of the ρ and a 1 mesons. The presence of the mass jump in fact depends non-trivially on the oscillation of the meson polarization function. This implies that the solution for the pole of the propagator will exist depends on the external parameters. For example, there can be certain values of the magnetic fields for which no mass jump may will occur (see for example Fig. 10(a) for eB = 0.05 GeV 2 ) within a certain range of temperature. Moreover, one should keep in mind that the polarization function also requires a regularization prescription. In our two step regularization scheme, the dimensional regularization is the essential first step to obtain the Lorentz structure for the vector and axial-vector mesons. As mentioned earlier, the Lorentz structure can not be achieved systematically in thermo-magnetic case with the cut-off procedure commonly used. Thus, it is very interesting to study the similar analysis in other covariant regularization prescription such as Paul-Villars method to conclude about the regularization scheme independent qualitative properties of the mesons . dx g µν ∆ + x(1 − x)q µ q ν (−2 ln z) + g µν ∆ + x(1 − x)q µ q ν (2ψ(z) + 1/z) (E36)