Dynamical Mass Generation in Pseudo Quantum Electrodynamics with Gross-Neveu Interaction at finite temperature

We study the dynamical mass generation in Pseudo Quantum Electrodynamics (PQED) coupled to the Gross-Neveu (GN) interaction, in (2+1) dimensions, at both zero and finite temperatures. We start with a gapless model and show that, under particular conditions, a dynamically generated mass emerges. In order to do so, we use a truncated Schwinger-Dyson equation, at the large-N approximation, in the imaginary-time formalism. In the instantaneous-exchange approximation (the static regime), we obtain two critical parameters, namely, the critical number of fermions $N_c(T)$ and the critical coupling constant $\alpha_c(T)$ as a function of temperature and of the cutoff $\Lambda$, which must be provided by experiments. In the dynamical regime, we find an analytical solution for the mass function $\Sigma(p,T)$ as well as a zero-external momentum solution for $p=0$. We compare our analytical results with numerical tests and a good agreement is found.


I. INTRODUCTION
In the last decades, the interest in studying twodimensional theories has been increased since the experimental realization of graphene [1,2] and other related materials, such as the transition metal dichalcogenide monolayers [3]. The main goal is to derive either new quantum phases of matter [4][5][6][7] or calculate renormalized parameters that change the electronic properties of these materials [8][9][10][11][12][13], opening possibilities for future technological applications (in particular spintronics [14], valleytronics [15], and electric-field tunning of energy bands [16,17]). Before this context, in the realm of high-energy physics, several works discussed the possibility of dynamical mass generation for massless Dirac particles, yielding a phase transition to a new quantum state of matter in which the chiral symmetry is broken [18][19][20][21]. This is generated due to the electronic interactions in the plane and it may occur even at finite temperatures. Because electrons in these materials obey a Dirac-like equation, therefore, some authors have also considered the realization of dynamical mass generation in these two-dimensional systems [22][23][24][25].
In this case, the fermionic field exhibits a massive phase due to the spontaneous symmetry breaking, described by the so-called gap equation. This is usually calculated at the large-N expansion, where N is the number of fermion fields and the gap equation is calculated at order of 1/N [29,32]. The dynamical mass generation is a typical effect of the nonperturbative regime, hence, it is common the application of the Schwinger-Dyson equations (SDEs) for calculating the critical parameters that describe the phase transition. The SDEs are an infinite set of coupled-integral equations, relating all of the interacting Green functions of the model [19,20,33]. Fortunately, within some approximative scheme, it is possible to find a truncated set of equations for calculating the desired Green functions, in particular, the electron-self energy that provides the chiral symmetry breaking. For pseudo-quantum electrodynamic (PQED) [23,34], also called reduced quantum electrodynamics [35], the dynamical mass generation has been studied both at finite temperatures [36] and at zero temperature with the presence of the Gross-Neveu interaction [31]. Nevertheless, the effect of the Gross-Neveu interaction in PQED at finite temperatures has not been considered until now.
In this work, we describe the dynamical mass generation in PQED coupled to a Gross-Neveu interaction at finite temperature. We use the Matsubara formalism in order to include the effect of the thermal bath into the Schwinger-Dyson equation for the electron. This is dependent on both the gauge-field propagator and the auxiliary-field propagator, obtained after we use a Hubbard-Stratonovich transformation into the four-fermion interaction. These two bosonic propagators are calculated in order of 1/N , which is consistent with our assumption of strong-coupling limit. Thereafter, we use this result into the Schwinger-Dyson equation for the electron self-energy and calculate the full electron propagator in the dominant order of 1/N in the nonperturbative limit. From this result, we conclude that a mass function is dynamically generated whether the number of fermions is less than a temperature-dependent critical parameter N c (T ). At large temperatures, we find that N c (T ) ≪ 1, hence, the dynamically generated mass vanishes and the system is in the gapless phase.
This paper is organized as follow. In Sec. II. we show our model and perform the large-N approximation. In Sec. III. we write the truncated set of Schwinger-Dyson equations within the unquenched-rainbow approximation. In Sec. IV. we calculate the mass function in the static regime and obtain the critical parameters for the phase transition. In Sec.V. we use the zero-mode approximation in the dynamical regime for calculating the mass function. In Sec.VI. we summarize and discuss our main results. We also include three appendices, where we give details about the angular integral, kernel expansion, and the numerical results for the mass function.

II. PSEUDO-QUANTUM ELECTRODYNAMICS WITH GROSS-NEVEU INTERACTION
We consider N fermion fields constrained to the plane, whose interaction is described by Pseudo-quantum electrodynamics (PQED) [34]. Furthermore, we assume that these particles involve a contact interaction, given by the Gross-Neveu (GN) action [37]. Therefore, in the Euclidean space-time, the action of the model reads where F µν is the field intensity tensor of A µ , which is the gauge field, ξ is the gauge-fixing parameter, ψ a is a four-component Dirac field, with the flavor index a = 1, . . . , N , the dimensionless coupling constant e is the electric charge, γ µ are the Dirac matrices in the four-rank representation, whose algebra is defined by {γ µ , γ ν } = −2δ µν , the coupling constant G 0 describes the strength of the GN interaction and has unit of mass in the natural system of units ( = c = 1). In the perturbative regime at zero temperature, the model in Eq. (1) for massive fermions has been used to describe the renormalization of the bandgap for WSe 2 and MoS 2 [13]. This result, nevertheless, requires a bare-mass term such as m 0ψa ψ a which is renormalized at oneloop and provides a beta function for the mass in terms of the RG scale. From this result, one concludes that the renormalized mass is dependent on the electronic density and a good agreement with experimental data has been found in Ref. [13] In the nonpertubative limit, the mass term is generated due to interactions even if we start with m 0 → 0 [31]. Here, we generalize this result by including a thermal bath of temperature T .
From Eq. (1), we find the gauge-field propagator, namely, and the fermion propagator Before we discuss the vertex interactions, let us apply the Hubbard-Stratonovich transformation, which converts the four-fermion interaction into a Yukawa-type interaction by including an auxiliary field ϕ and the coupling constant g = G 0 N . In this case, we replace the four-fermion interaction by the following scheme where L GN = −G 0 (ψ a ψ a ) 2 /2. Therefore, after that we find Equation (5) must be supplemented by the motion equation of ϕ, given by which proves that the transformation does not change the dynamics of the model at classical level. Furthermore, the bare auxiliary-field propagator is and Eq. (1) reads The Yukawa-type vertex function is easily obtained from Eq. (8) and is given by −1/ √ N . On the other hand, for summing the self-energies in the large-N expansion for PQED, we shall replace e 2 → λ/N , where λ is taken fixed at large N . This allow us to sum over all of the diagrammatic contributions in order of 1/N , which is an infinite sum, unlike the standard perturbation in e. Therefore, the PQED vertex reads λ/N γ µ , describing the electromagnetic interaction.

III. TRUNCATED SCHWINGER-DYSON EQUATION AT FINITE TEMPERATURES
In this section we present the Schwinger-Dyson equations, obtained from Eq. (8), that describes the quantum corrections for the two-point functions. In principle, this is a very complicated set of integral equations for all of the full Green functions of the model. From now on, we assume the ladder approximation, also called rainbow approach [38], which consists of neglecting quantum corrections to the vertex functions. It is worth to note that this may be corrected by the Ball-Chiu vertex [39] when one wishes to preserve the Ward-Takahashi identity.
Because we are interested in dynamical mass generation for the electrons, we need to find a closed solution for both gauge-field and auxiliary-field propagators. These shall be given by the large-N approximation, thus the bosonic-field propagators are calculated in the unquenched approximation [40]. The left-hand side is the inverse of the full propagator of the auxiliary field, the first term in the right-hand side corresponds to the bare auxiliary-field propagator, and the second term Π(p) is the quantum correction. The function Its analytical expression is given by [31] where ∆ −1 ϕ is the full propagator of ϕ and Π(p) is given by the fermionic loop, hence where S F (p) is the full fermion propagator. In the lowest order of 1/N , we find that Π(p) = − p 2 g 0 , with g 0 = 1/4 [21]. Therefore, the full propagator of the auxiliary field is given by In Fig. (2), we show the diagrammatic representation of the Schwinger-Dyson equation for the gauge-field propagator [31].
The Schwinger-Dyson equation for the gauge-field propagator. The left-hand side is the inverse of the full propagator of the gauge field, the first in the right-hand side corresponds to the bare gauge-field propagator, and the second term Πµν (p) is the exact vacuum polarization tensor. The function Γ µ [p, k] → γ µ corresponds to the approximated vertex.
Its analytical expression is given by where ∆ −1 µν is the full propagator and Π µν is the polarization tensor due to the electromagnetic interaction. This is given by Similarly to the previous case, in the lowest order of 1/N , we find Π µν = λ p 2 /8P µν , where P µν = δ µν − p µ p ν /p 2 . Therefore, the full propagator reads Using Eqs. (2) and (13) in Eq. (14), we find in the Landau gauge, i.e, with ξ = ∞.

C. Matter field
In Fig. (3), we show the diagrammatic representation of the Schwinger-Dyson equation for the fermions field [31]. Its analytical expression is given by where S F is the full fermion propagator and Ξ is the electron self-energy, given by where Ξ λ is the PQED contribution and Ξ g is the GN contribution, namely, In order to find the full fermion propagator, it is convenient to decompose this Green function into its irreducible parts. Hence, we adopt the following ansatz where A(p) yields the wavefunction renormalization and Σ(p) is called mass function. Note that Σ(p) = 0 implies a dynamical mass generation [21]. In order to obtain the mass function, we take the trace operation in both sides of Eq. (16). Thereafter, we substitute Eqs. (3), (18), (19) and (20), within the rainbow approximation, On the other hand, for calculating the wavefuntion renormalization, we multiply Eq. (16) by / p and, after calculating the trace operation, we obtain It is clear that Eqs. (21) and (22) are a coupled set of equations for Σ(p) and A(p). However, in the dominant order 1/N , we may replace A(p) → 1 in Eq. (21). This is also in agreement with the rainbow approximation.
Next, we use the Matsubara formalism for introducing temperature effects into Eq. (21). In the imaginarytime formalism, the main step is to change the timecomponent integrals into a sum over Matsubara frequencies [41], i.e, we must perform with and ω n = where From now on, we assume that Σ m (p, T ) = Σ(p, T ) and, therefore, the mass function only depends on the dominant vibrational mode m = 0. The standard procedure for calculating the mass function is to convert the integral equation into a differential equation for Σ(p, T ). It turns out that there exist two main operations, namely, the sum over n and the angular integral that must be performed before finding this result. We shall explore some approximations for these operations.
Before performing further approximations, let us show that the sum over the vibrational modes n is convergent. The convergence of this sum for PQED at finite temperature has been shown in [36]. In Eq. (25), the sum over n may be written as where we conclude that (A 2 , C 2 ) > 0 and, because we consider g > 0 [37,42], we also find that B > 0. Next, we apply the Cauchy integral test [43], where the sum is written as with u (+) n = (u n +u −n )/2. The test imposes that whether u(x) = u x is positive, continuous, and decreasing in the interval [1, ∞], hence, the integral over u(x) is finite and, therefore, convergence is derived. This is exactly our case and one may easily conclude that Therefore, Eq. (30) is a convergent series for all n.
IV. STATIC REGIME p0 = 0 The static regime, also called instantaneous-exchange approximation [44], consists of neglecting all of the timecomponents of the bosonic-field propagator in Eq. (21), with the consideration that the interaction vertex being only γ 0 in Eq. (18). This is usually realistic in cases where the electron velocity is much less than the light velocity. In these systems, it is possible to show that charged particles interacts through the Coulomb potential. This is given by the Fourier transform of the gauge-field propagator in Eq. (2) after we perform p µ = (p 0 , p) → (0, p), where we will use the following notation (0, p) ≡ (0, P ). Next, let us consider both zero and finite temperature cases for calculating Σ(P, T ).

A. Zero Temperature Case
After considering the static regime with T → 0 in Eq. (25), we find where and ∆(P − K) = 1 The integral over k 0 in Eq. (34) is easily solved and the angular integral may be solved by using the identity where K(x 0 ) is the complete elliptic integral of the first kind with x 0 ≡ −4P K/(P − K) 2 . Note that the angular integral for the gauge-field term is very similar and that we will neglect the terms in order of O[(gg 0 ) −2 ]. Replacing Eq. (37) in Eq. (34), we find where the kernel K(K, P ) reads Next, we approximate the kernel for its lowest order terms, namely, Note that Θ(P − K) is the standard step function. Using Eq. (40) in Eq. (38), we have On the other hand, by taking derivatives in respect to P in both sides of Eq. (41), we obtain a differential equation for the mass function, given by where the critical number of fermions is given by Eq. (42) is also supplemented by two equations, namely, lim P →0 which provides the infrared behavior and lim P →Λ which shows that the mass function is expected to vanish at large P . Using the approximation P/ P 2 + Σ 2 (P ) ≈ 1 in Eq. (42), which holds for large-external momentum, we find the analytical solution where γ = −1/2 + i/2 N c /N − 1 and (A 1 , A 2 ) are arbitrary constants. In the dynamical regime, it has been show that the dynamical mass generation only occurs for N < N c [31]. The same conclusion form Eq. (43) is also obtained from the static regime. Using Eq. (43) in γ, λ = 4πN α c , and solving N c (α c )/N = 1 for α c , we find the critical fine-structure constant, given by in terms of fine-structure constant the mass generation condition is α > α c .

B. Finite Temperature
In order to describe temperature effects, we include the Matsubara frequencies in Eq. (34). Thereafter, we solve the angular integral using Eq. (37) while for the Matsubara sum we use the identity where ǫ K ≡ K 2 + Σ 2 (K) and is the Fermi-Dirac distribution. After using these properties, we find where the temperature-dependent kernel K β (K, P ) reads Next, following the same steps as in the previous section, we find the differential equation for the mass function, namely, The linearized version of Eq. (52) may be obtained at large-external momentum and by doing P → Λ in h β (P ). This allow us to obtain a temperature-dependent critical number N c (T ), given by where N c (0) ≡ N c = (C 1 /2 + C 2 )/π is given by Eq. (43) in terms of λ. On the other hand, for large-T , we obtain N c (T ) → 0, which implies that there exist no dynamical mass generation for any λ.
In the large-external momentum, the solution of Eq. (52) is given by Eq. (46), after we replace N c → N c (T ). Hence, we obtain a temperature-dependent exponent γ(T ) = −1/2 + i/2 N c (T )/N − 1. Therefore, we find where the arbitrary constants (A 1 , A 2 ) may be dependent on the temperature. At the critical point N = N c (T ), we have γ(T ) = −1/2. This critical point also may be described in terms of α, using Eq. (43) and that e 2 = 4πα = λN . Indeed, after a simple calculation, we obtain a temperature-dependent critical coupling constant α c (T ), namely, Because the ratio N c (T )/N is a monotonically increasing function of α, it follows that the mass generation only occurs for α > α c (T ). In Fig. 4, we plot both N c (T ) and α c (T ) and show that whether the temperature increases, hence, the system quickly goes to a gapless phase. Note that for small enough temperatures, the critical coupling constants remain almost unchanged.
We also may find numerical results (see App. A for more details) of the integral equation for Σ(P, T ) in Eq. (50). In Fig. 5, we plot these results and compare this with our analytical solution given by Eq. (55). V. DYNAMICAL REGIME p0 = 0 In this section we recover the retardation effects by assuming p 0 = 0 in Eq. (21). For T = 0, this regime has been discussed in several contexts.
After a simple changing of variables, the angular inte-gral of the PQED sector (see details in App. B) reads where, and Similarly, for the GN sector, we find where K(x n ) and Π(y n , x n ) are the elliptic integral of the first and third kind, respectively. In particular, note that for g = 0, we have that Π(0, x n ) = K(x n ), hence, the kernel of the GN sector vanishes, as expected. The integral representation for these functions are and Π(y n , . (63) The resulting kernel, after using the angular integral, is a very complicated function of n. This preclude us of solving the whole Matsubara sum, as we have did in the static regime. Therefore, we shall investigate the lowest order contribution, provided by the zero-mode n = 0 approximation.

VI. ZERO-MODE APPROXIMATION
In this section we consider the fundamental vibrational mode n = 0 in Eq. (25). Using Eqs. (62) and (63), we find For calculating an analytical solution of Eq. (64), we write the elliptic integral as the hypergeometrics Appell F 1 . Because F 1 has an expansion in its parameters, we may find a simplified version of the kernel (see App. C for more details). Hence, Eq. (64) yields where f (k) ≡ kΣ(k)/ (πT ) 2 + k 2 + Σ 2 (k) . By taking derivatives in respect to p in both sides of Eq. (65), we find with a = (C 1 +C 2 )T /2π and b = C 2 T /2π, both constants with dimension of mass. Next, we neglect the nonlinear terms by using (πT ) 2 + p 2 + Σ 2 (p) ≈ (πT ) 2 + p 2 , which is expected to be reasonable at large p [36,45]. After deriving Eq. (66) in respect to p, we obtain the differential equation which is supplemented by two conditions for the infrared and ultraviolet regimes. These are given by lim p→Λ Σ(p) = 0 (68) and lim p→0 p 2 dΣ(p) dp + bΣ(p) (πT ) 2 + p 2 = 0.
Unfortunately, Eq. (67) has not an analytical solution for arbitrary values of the set (p, T ) of variables. Nevertheless, we may consider a possible linearized version for large external momentum.

A. Large-external-momentum approximation
In this section we investigate the linearized version of Eq. (67) by assuming p ≫ T, Σ(p). In this case, the behavior of the generated mass is given by which yields where 1 F 1 is the confluent hypergeometric function and G is the Meijer G-function [46]. The sign of this Meijer G-function changes for different values of the external momentum. Because we have not observed such behavior in the numerical results of Eq. (64), from now on, we shall take A 2 = 0 for the sake of agreement with the integral equation. In Fig. (6), we plot our analytical solution in Eq. (71) and compare this with the numerical results of the integral equation in Eq. (64).

B. Zero-external-momentum approximation
In this section, we consider the limit Σ(p = 0, T ) = m(T ) in Eq. (25), using the zero-frequency mode. In this case, the angular integral provides a constant factor of 2π. Hence, we find After calculating the integral over k, we have Here, nevertheless, we define that C 1 T and C 2 T are fixed in order to be in agreement with the fact that m = 0 for T = 0. In Fig. 7, we plot the numerical solutions of m(T ) for both g = 0 and g = 1.0, which proves that the GN interaction increases the gapped phase. Furthermore, in Fig. 8, we show that the GN interaction increases the critical temperature, as expected.

VII. SUMMARY AND OUTLOOK
In this work, we investigate the influence of temperature on the dynamical mass generation for fermions in PQED with a Gross-Neveu interaction. In order to do so, we used the Schwinger-Dyson equation for the electron in the dominant order 1/N , negleting quantum corrections to the vertices at the finite temperature. In the static regime, we were able to solve the whole Matsubara sum. Thereafter, we calculated the critical coupling constant α c and the critical numbers of fermions N c as a function of T and the cutoff Λ. From these results, we concluded that the system goes to a gapless phase whether we increase the ratio T /Λ. Furthermore, we obtained an analytical solution for the mass function Σ(p), which is in good agreement with our numerical results.
In the dynamical regime, we used the zero-mode approximation in the sum over the Matsubara frequencies, which allowed us to calculate the analytical solution Σ(p). This function agrees with our numerical results for largeexteral momentum, while for small-external momentum some deviations have been found. Moreover, we found a gap equation that provides the value of the mass function in the zero-external momentum approximation, i.e, Σ(p = 0, T ) → m(T ). In this case, the numerical tests show that the critical temperature increases as we increase the strength of the fermionic self-interaction.
Several generalizations of this work may be perfomed. For instance, one would investigate the main effects of vertices quantum corrections into the critical parameters for the mass generation. In particular, it would be relevant to find a way to increase the gapped phase for finite temperature. In this phase, the competition between α and g yields a tunneable mass which may be relevant for applications in two-dimensional materials [13]. In this Appendix we briefly review the main steps in order to find our numerical results. For the sake of simplicity, note that Eq. (25) is a kind of Fredholm integral equation of first kind, defined by where K(x, y, Z(y)) is the kernel of the integral equation and Z(x) is an unknown function. Throughout this work we convert this into a differential equation, allowing us to find analytical solutions for each approximation, namely, Eq. (55) in the static approximation and Eq. (71) in the zero-mode approximation. These solutions, nevertheless, are not expected to well describe the behavior of the mass function for small external momentum.
For properly describing the behavior of Σ(p) for any p, we may find numerical results of both Eq. (50) and Eq. (65). It is easy to conclude that these equations are essentially a kind of Fredholm integral equation whether we adjust the parameters, variables, and kernel. Next, let us summarize our main steps. First, we use the trapezoidal quadrature rule for calculating the integral over the kernel, hence, b a dyK(x, y, Z)Z(y) → h 2 where h is the size of each interval and y i is the iterative variable (plays the role of loop momentum k).