High-temperature expansion of the grand thermodynamic potential for scalar particles in crossed electromagnetic fields

The problem of a scalar particle in a constant crossed electromagnetic field ($\mathbf{E}\perp\mathbf{H}$ and $|\mathbf{E}|=|\mathbf{H}|$) is examined. The high-temperature expansion of the grand thermodynamic potential and vacuum energy with account for non-perturbative corrections are derived. The contribution from particles and antiparticles is considered separately. It is shown that the non-perturbative corrections depend on boundary conditions but do not depend on fields.


Introduction
The present work is a continuation of the study devoted to the high-temperature expansions of the grand thermodynamic potential [1][2][3][4]. Thus, in the paper [3] we derived the high-temperature expansion of the grand thermodynamic potential in the case of a constant homogeneous magnetic field (F F = 2B 2 ,F F = 0) as well as the explicit form of the non-perturbative corrections. As a next step we see the analysis of a constant crossed electromagnetic field, when E and H vectors are orthogonal and equal in absolute value. There are several reasons why this field configuration is interesting. Firstly, all fields appear to be crossed for an ultrarelativistic particle in the co-moving frame of reference. Secondly, both electromagnetic invariants are trivial, F µν F µν =F µν F µν = 0.
One of the main results of the series of papers [1][2][3][4] is the general asymptotic formula for the expansion of the Ω−potential (see details of derivation in [4]). As for particles obeying Bose-Einstein statistics, The ≃ sign indicates that the expansion in ascending powers of β is asymptotic, and the terms exponentially suppressed in temperature (β → 0) are discarded. Space-time dimension D := 4, and the term with l = −1 should be considered as a limit, i.e. ζ(1)/Γ(0) = −1. It is seen that at some values of k and n the first terms may possess singularities as ν → 0 (coefficients ζ + k (ν) are always regular), which must be canceled exactly by singularities coming from σ l ν (µ). The functions σ l ν (µ) and the coefficients ζ + k (ν) entering the expansion (1) are determined by zeta function constructed by means of Laplacian type operator H(ω) where C-contour runs upwards a little to the left of the imaginary axis. The operator H(ω) is a Fourier image over time of the Klein-Gordon type operator and possesses a spectrum bounded from above. The "+" index of zeta function reminds us that its values are determined only by positive eigenvalues of H(ω). The coefficients ζ + k (ν) are the coefficients of the asymptotic expansion of zeta function for large ω and the functions σ l ν (µ) are determined in the following way σ l ν (µ) = ∞ 0 dω(ω − µ) l ζ + (ν, ω).
It is the functions σ l ν (µ) that contain exponentially suppressed in fields corrections and which calculation is the most difficult.
The contribution of antiparticles to the thermodynamic potential is derived from (1) by changing the sign of a chemical potential and by a simultaneous replacement ζ + k (ν) → ζ − k (ν), where the coefficients ζ − k (ν) are determined from the expansion It should be noted here that for the configuration of fields and plates under study the coefficients ζ + k (ν) coincide with ζ − k (ν). Hereinafter both types will be denoted by ζ k (ν). For the field configuration at issue the problem is complicated by the fact that the expression for the heat kernel (exponent on the right-hand side of (2)) derived in [5] (see also [16]) is not appropriable. The standard definition of the vacuum state in a stationary field can be used only when the work of the electromagnetic field is smaller than the energy required to create a particle-antiparticle pair from the vacuum, ∆A 0 < 2mc 2 , see [6] (the situation with strong enough electric fields is considered in, e.g., [7]). Therefore, we are forced to place constraints on a space extension of a system, with the result that a problem acquires nontrivial boundary conditions. In this respect, the existence and the form of non-perturbative corrections are of major interest.
Section 2 is devoted to the formulation of the problem, the computation of the spectral density, and the derivation of the valuable relations between the parameters of the theory that will be used when calculating σ l ν functions, which define the one-loop contribution to the grand thermodynamic potential. In Section 3, the explicit expressions for the first six coefficients ζ k (ν) are found. The major trick used there is that zeta function can be represented as an integral of the function defining the spectrum of the problem. The alternative method of calculation of the coefficients is presented in Appendix A. Section 4 is devoted to the calculation of σ l ν functions. Despite the fact that an explicit calculation of the functions at arbitrary l is impossible, emerging structures allow one to perform an exact computation in the case of non-negative integer l. The method described allows to derive the expressions for any l, explicit calculations were carried out up to l = 3. The last section presents the explicit expressions for finite and divergent parts of the grand thermodynamic potential and the renormalized vacuum energy taking into account the non-perturbative corrections. The non-perturbative corrections turn out to be independent of the electromagnetic fields and are exponentially suppressed at large mL, where L is the extension of the system along the field E.

Spectrum
Let us consider the eigenvalue problem for the Klein-Gordon operator with the constant homogeneous crossed electromagnetic field where A µ is given in the Coulomb gauge and, for definiteness, E > 0. The naive definition of a particle is applicable in this field provided EL < 2m, where L is the size of the system along the z axis (see for details [6]). Therefore, we impose the zero Dirichlet boundary conditions on the wave function and consider the problem on the segment z ∈ [−L/2, L/2]. Separating the variables, we obtain where the particle charge is included into the definition of the electromagnetic potential. This differential equation is reduced to the Airy equation and has the general solution The spectrum ε is found as the solution to the equation Obviously, the spectral density with respect to ε is equal to Further, we shall need the inequalities following from (7). Averaging (7) with respect to the eigenstate, we find for particles (ω > 0) Then, for ε ≥ 0, Using the first inequality, we deduce from (7) that The spectral density is zero where the above inequalities are not satisfied for ω > 0.

Coefficients ζ k (ν)
Despite the fact that the spectral equation (9) cannot be solved explicitly, it is possible to represent the spectral zeta function as an integral of a function defining the spectrum (the so-called Gelfand-Yaglom formalism, see, e.g. [8]) where the contour γ runs along the imaginary axis downwards, and p y in h − , h + is set to zero. Let us prove that ζ + k (ν) = ζ − k (ν). It is easy to show that the spectral equation for antiparticles coincides with (9) with the replacement p x → −p x . Therefore, ρ(ε; −ω, p x , p y ) = ρ(ε; ω, −p x , p y ). Due to this relation To calculate the coefficients ζ k (ν), it is necessary to expand zeta function into a series for large ω, which corresponds to large negative h ± . Using the asymptotic expansion of Airy functions [9], we get Then, for the logarithm of the expression (18), we have It is taken into account that the cotangent tends to ±i on the upper and lower parts of the contour. The last two logarithms in (20) are connected by the replacement h + ↔ h − . It is easy obtain the first several terms of the expansion of the logarithm: As we will see later, the terms presented suffice to calculate ζ k up to k = 6. The fact that we are interested in the explicit form of the first six coefficients is attributable to that we managed to calculate the functions σ l ν till l = 3. Thus, we shall know the expansion of the Ω-potential up to β 3 (see the formula (109)). It should be noted that in order to derive the finite and divergent, as β → 0, parts of the expansion, it is sufficient to know the coefficients till ζ 4 . Consider the first contribution Enclose the contour on the cut of the function ε 1/2−ν and perform the integral over ε: Here we have introduced the notation ω ± := ω ± EL/2. Note that the contributions with E and −E cannot be considered separately as the integral over ε is divergent at any value of ν. After integrating over p x , we arrive at The hypergeometric function should be expanded in the vicinity of unity. The easiest way to do this is by using the relation It is not difficult to see that only the second term in (25) gives the leading contribution as ω → +∞, so Expanding the derived expression into a series in 1/ω, we get , Consider the term The contributions with E and −E can be treated separately. It is convenient to enclose the contour γ to the right and calculate the integral using residues. After being integrated, the contribution to ζ + (ν, ω) reads as Then, the contribution to ζ k from (28) becomes It should be pointed out that the expressions (26) and (29) exactly coincide with the known answer for zeta function as E → 0: In the contribution from the last line of (21) we integrate by parts It is sufficient to consider only the integrals of the form is the binomial coefficient, and the sum over n is finite as − 2 3 α is a positive integer number. As a result, the contribution to ζ + (ν, ω) reads Using this formula for the third term in (21), we obtain Expanding this expression into a series in 1/ω, we arrive at For the fourth term in (21), we get the following contribution to zeta function The contribution to ζ 5 Collecting all the contributions together and taking the coefficients into account, we find The above mentioned expressions for ζ k (ν) allow one to conjecture the general structure at arbitrary k: The rigorous proof of (44) follows from the analysis of the expressions (26),(29), and (36). An alternative method of calculation of the coefficients by use of Dyson series is represented in Appendix A.

Functions σ l ν (µ)
In accordance with the general formulas (see Introduction), for the high-temperature expansion to be obtained, one needs to find the expression for the function where S := L x L y , in the form of an analytic function of ν and l. It was shown in [1] that the integral over ε converges when Re ν < 1.
The integral over ω converges for Therefore, it is useful to calculate (45) in the region (46), (47), where the multiple integral (45) converges, and then to continue σ l ν (µ) by analyticity to the required "physical" values of the parameters ν and l. Recall that, according to the Hartogs theorem (see, e.g., [10]), the function that is analytic with respect to each variable is analytic with respect to all of them. The uniqueness of analytical continuation also holds for such functions. As follows from the general analysis given in [1], the function σ l ν (µ) is a meromorphic function of ν and l. It possesses the singularities in the form of simple poles at provided that there is a neighborhood of the point ω = 0 that does not contain the points of the particle's energy spectrum and µ belongs to this neighborhood. First, we integrate (45) over p y . After the replacement ε → ε − p 2 y , the spectral density ρ(ε) becomes independent of p y , and the integral over p y is reduced to Further, we make the integration variables dimensionless introduce convenient notation and pass to the light-cone variables The region of integration with respect to these variables is determined by the inequalities (12).
we arrive at As seen from this expression, we can integrate over the variable v as it was done above with the variable p y .
To this aim, we shift the integration variable Then the integrand of (54) includes On making the redefinitions (50), (52), (53), and (55), the inequality (13) has the form Therefore, the first θ-function in (56) can be removed. As a result, the integral over v becomes and taking into account the inequality (57), we obtain Further simplification of this integral is impossible without knowledge of the explicit expression for the function ϕ. It appears at first sight that (60) possesses singularities at which contradicts the general statement (48). However, for such values of ν and l, the integral on the second line of (60) understood in the sense of analytic continuation goes to zero (see (90)) and, consequently, σ l ν (µ) 4.1 Functions σ l ν (µ) for nonnegative integer l The presence of Γ(−l) in the denominator of (60) allows one to obtain the exact expression for σ l ν (µ) at l = 0, ∞. To this end, it is sufficient to investigate the singularities of the integral on the second line in (60) in the complex l plane near l = 0, ∞. As the theorem 1 shows, these singularities are the poles with the residues found by expansion of the integrand in the asymptotic series. There is no need to evaluate the integral. Namely, (see for details [11][12][13]) Theorem 1. Let ϕ(x) be absolutely integrable on (0, Λ] and the following asymptotic expansion takes place for x → +0. Then the function is analytic for Re λ > −1 and can be analytically continued to the region Re λ > −2 − N , where it possesses the simple poles at the points λ = −k, k = 1, N + 1, with the residues a k−1 , respectively. In order to apply the theorem to (60), we pass from the integration variable ε to where it is supposed that µ + w < 1. Let us stretch the integration variable Then the integral on the second line of (60) can be cast into the form where and The last expression is obtained from the solution of the equation h − = 0. The function g(ξ) is bounded on the integration interval except possibly the vicinity of the point ξ = 0. Therefore, in order to evaluate (60) at l = 0, ∞, it is sufficient to derive the asymptotic expansion of the integrand for ξ → +0.
For ξ → 0, it is useful to split the integration region of the variable u in (67) into the three intervals: where The boundaries of the integration region (71) are obtained from the solution of the equation ε = ε 0 , where ε 0 is some constant independent of ξ. For the intervals A and C, the quantity ε ∈ [ε 0 , 1]. For the interval B, both h + and h − are nonnegative as long as ε < ε 0 .
Intervals A and C. On the interval A, the integration variable u → c −1 for ξ → 0. In this case, Eq. (9) has the form where ε ∈ [ε 0 , 1]. Then it follows from (72) that However, u ≈ c −1 . Consequently, for ξ → 0, Eq. (9) does not possess solutions on the interval A.
On the interval C, the integration variable u ∼ ξ −2 for ξ → 0. Hence, it is useful to redefine the integration variable u → ξ −2 u. Then the integration limits become for ξ → 0. On stretching u, the additional factor ξ 2ν appears in the expression (67) for the function g(ξ) (on the interval C). If Eq. (9) has solutions on the interval C, then the integration over u in (67) is removed and u = 4w 2 /c + o(1), where o(1) does not contain the powers of ξ ν . Then, in developing g(ξ) as a series in ξ, the factor ξ 2ν cannot be canceled out, i.e., the expansion of g(ξ) (on the interval C) in the vicinity of ξ = 0 does not contain integer powers of ξ. Consequently, as follows from the theorem 1 and the form of the integral over ξ in (67), the interval C does contribute to the poles at nonnegative integer powers of l.
Interval B. On the interval B, the arguments of the Airy functions entering into the equation specifying the spectrum (9) are negative and tend to −∞ for ξ → +0. Therefore, we can employ the asymptotic expansion of the Airy functions for large negative arguments and solve Eq. (9) with respect to u, bearing in mind that ξ → 0. The first three terms of the expansion with respect to ξ are written as where N (ξ) → ∞ for ξ → +0. For example, in the leading order, we obtain from (9) that The solution of this equation reproduces the leading term of the expansion (75). The integral over u in (67) is reduced to the sum over the roots of Eq. (9) since on the interval C and for ξ → +0: Substituting the expression (75) for u n into the resulting sum and expanding the outcome in a series with respect to ξ, we have Thus, employing the theorem 1, we deduce The explicit expressions for the sums over n are presented in Appendix B. Using these expressions, we can write The functions T k α (4c) are exponentially suppressed for large c. We see from formulas (45), (81) that σ 1 ν (0) vanishes. In general, for the system at issue Indeed, it follows from (68) that for µ = 0. As long as Eq. (9) remains unchanged under (83) and ϕ(h + , h − ) = ϕ(h − , h + ) (see (77)), the expansion of g(ξ) with respect to ξ contains only even powers of ξ. Therefore, we have (82). It is clear from (4) that Then formula (82) implies that σ 2k+1 Also we obtain Taking into account the inequalities (14), it is not difficult to check that the contribution of antiparticles to σ l ν (µ) is equal to the contribution of particles with µ → −µ. Hence, the high-temperature expansion of the thermodynamic potential including the contributions of particles and antiparticles does not contain σ 2k+1 ν (µ). The terms with σ 2k ν (µ), k = 1, ∞, do not contribute to the high-temperature expansion either (see (1)). Consequently, the one-loop contribution to the Ω-potential coming from particles and antiparticles does not include the second term with l = 1, ∞. Notice that, in contrast to the system of charged particles in a homogeneous magnetic field (see [3]), the functions σ l ν (µ), l = 0, ∞, do not contain the terms that are nonanalytic with respect to the coupling constant or the external field. Such contributions do not arise in perturbative in ξ solution of Eq. (9) with respect to u. As for σ l ν (µ), l = 0, 3, this property is seen directly from the expressions (81), (86).

Function σ −1 ν (µ)
In order to find the high-temperature expansion for the thermodynamic potential of bosons, we need to derive σ −1 ν (µ). The considerations of the previous subsection are not applicable in this case. We did not succeed in finding the exact expression for σ −1 ν (µ), but we managed to find it under the assumption that (see (51)) Let us change the integration variable in the integral on the second line of (60), and rewrite it as a contour integral where the contour C goes from −∞ a little bit lower than the real axis, encircles the origin, and then runs to −∞ a little bit higher than the real axis, h ± = c 2/3 u 1/3 (ε ± 1), and the principal branches of the multivalued functions are taken. The logarithmic derivative entering into (89) does not possess singularities out of the negative part of the real axis and tends to zero for |ε| → ∞. Therefore, taking Re ν sufficiently large, we can deform the contour C and reduce the integral (89) to the integral over the branch cut of the power function in (89). As a result, we obtain where Bearing in mind the conditions (87), it is not difficult to see from the asymptotic behavior of the Airy functions that, in the given integration region, is exponentially small. Therefore, up to exponentially suppressed terms at (c/w) → ∞, we have Unfortunately, even in this case, we did not succeed in exact evaluation of the integral (90). Below, we shall derive the leading in c/w contribution to σ −1 ν (µ), the terms diverging for ν → 0, and take into account the asymptotics of σ −1 ν (µ) in the vicinity of the singular point in the complex µ plane where the chemical potential approaches to the lowest energy level. The integral (90) possesses the branch point in the µ plane when the pole of the logarithmic derivative tends to the point ε = 0.
We start with the leading contribution and the terms singular at ν → 0. Developing the logarithmic derivative as an asymptotic series, we deduce whereh ± := c −2/3 u −1/3 h ± . This asymptotic expansion has the same form for both the exact logarithmic derivative and the approximate expression (93). The terms presented in (94) are sufficient to find all the contributions diverging at ν → 0 for l = −1. On substituting the expansion (94) into (90), the integrals over ε of every term in the series are reduced to the beta-function: Having stretched the integration variable u → 2wu, the integral over u becomes (see [14]) Then the leading term of the expansion (94) gives the contribution to σ −1 ν (µ) of the form whereμ :=μ + Ez/m and, for convenience, the expression is written as the integral over z. As for the next term of the asymptotic expansion (94), we will obtain only the contribution that is singular at ν → 0. Using the above integrals, it is easy to see that this contribution to σ −1 ν (µ) is the pole part of Thus, it only remains for us to find the asymptotics of σ −1 ν (µ) in the neighborhood of the singular point of the µ plane. To this aim, we extract the pole contributions from the logarithmic derivative (94) that are closest to ε = 0: where Ai(−r − ) = Bi(−r + ) = 0 and r − ≈ 2.34, r + ≈ 1.17.
Substituting (99) into (90), integrating over ε, and stretching the integration variable u → 2wu, we obtain for the contributions of each of the poles. The integral (101) is singular when the expression in the square brackets vanishes for some u = u 0 ≥ 0. The main contribution to this singularity comes from the vicinity of the point u = u 0 . Expanding the expression in the square brackets near this point and keeping only the leading terms, we have where a ± := r ∓ (2w/c) 2/3 . The integral (101) diverges when s ± :=μ ± w approach the point Substituting the expansion (102) into the integral (101) for l = −1 and integrating, we come to in the leading order in (w/c). Keeping only the singular term at s ± → s 0 (a ± ) and replacing s 0 (a ± ) → 1, which is justified in the leading order in (w/c), we deduce where, in the last two terms, it was taken into account that ν → 0. The pole at ν → 0 in the expression (104) was discarded as its contribution is taken into account in the last term in (105). The finite terms at ν → 0 in (104), which are not singular for s → s 0 , can also be disregarded to the accuracy of the approximations made. For w/c ≪ 1, the expression (105) gives a good approximation for σ −1 ν (µ) and its derivatives with respect to the chemical potential on the interval |μ| < 1. The expression (105) has a branch point atμ = 1 and becomes complex forμ > 1. However, the exact expression for σ −1 (µ) is real-valued and well-defined for s ± < s 0 (a ± ). That is, there exists the interval of values of the chemical potential, which shrinks to zero for w/c → 0, 1 <μ + w < s 0 (a + ), where the expression (105) is not applicable. One can improve the expansion (94) such that, having integrated over ε and u, the leading contribution will approximate σ −1 ν (µ) uniformly on the whole interval of physically acceptable values of µ. To this end, one needs to put in the arguments of the Airy functions and employ the asymptotic expansion in (94) supposing that the expression in the parenthesis is large. Then, on evaluating the integral over ε, the integral (96) is replaced by which is just the Mellin transform. In the present paper we will not investigate the expression resulting from this procedure. Let us point out some other properties of (105). The term with logarithmic singularity in (105) gives the leading contribution to the total charge and the average number of particles at s ± → 1 in spite of the fact that this contribution is suppressed by the factor w/c in comparison with the first term. The first term in (105) can be obtained from σ −1 ν (µ) for a free particle [15]. One just needs to replaceμ →μ and integrate the resulting expression over z (see [16]). As we have already noted above, the contribution of antiparticles is obtained by replacement µ → −µ in the corresponding expression for particles. The contribution of antiparticles to the thermodynamic potential cancels out all the terms in (105) that are odd with respect toμ.
5 High-temperature expansion

Cancellation of singularities
The expression in square brackets in formula (1) is an entire function of the parameter ν, see [2]. This fact can be exploited for the indirect verification of the expressions for σ l ν whose explicit form is known for l = 0, 3. In the expansion (1) the terms with singularities at ν → 0 are It has been taken into account that there are no singularities in the first and in the second terms at β 2 and also that ζ 5,7 (0) = 0, see (44). It is easy to convince oneself that all the singularities are canceled in this case. This serves as an indirect verification of the fact that the divergent at ν → 0 (volume) terms in σ l ν (µ) were calculated correctly, e.g., the term with Γ(ν − 2) in σ 0 ν . Let us extract the terms with singularities at ν → 1/2: It has been taken into account that ζ 4,6 (1/2) = 0, see (44). It is not difficult to check that all the singularities are canceled as well. This suggests that the divergent at ν → 1/2 (surface) terms in σ l ν (µ) were calculated correctly, e.g., the term with Γ(ν − 3/2) in σ 0 ν .

Grand thermodynamic potential
Now we have everything to write down the explicit expression for the high-temperature expansion of the one-loop Ω-potential. Formulas (43), (81), (86), and (105) suffice to obtain the expansion up to β 3 terms.
We present here only the expression for the divergent and finite terms as β → 0. Thus, for bosons If we take into account the contribution from antiparticles, then The expression for the terms at positive powers of β of the high-temperature expansion (of the one-loop Ω-potential) with the contribution from antiparticles reads as It is easy to see that the derived expression is regular in the ν-plane; if there are singularities in gamma functions, they are annihilated by zeros of zeta function or the coefficients ζ k (ν), see (43). The vacuum contribution is not taken into account in (111), (112), and (113). The derived formulas can be applied to the systems in a local thermodynamical equilibrium. Notably, let the system be represented as a collection of subsystems with characteristic size L, which are large enough for the statistic description to apply, and on a scale L the electromagnetic field can be approximated as crossed, constant and homogeneous in a co-moving frame of reference. It is also assumed that the subsystems possess small acceleration in the co-moving frame of reference so as there is no impact of the metric on thermodynamic properties of the subsystem. Then, for such a system, we find from (111) that in the leading order of w/c. Here, T := β −1 , the chemical potential is incorporated into the definition parameter L entering (114) and, consequently, the entire coefficient at the term in question have only an order-of-magnitude estimate. Taking into account the contribution from antiparticles, we deduce (115)

Vacuum energy
In the paper [4] the high-temperature expansion of the grand thermodynamic potential for particles obeying the Fermi-Dirac statistics was derived Here, in contrast to the formula (1) There is an interesting fact that the high-temperature expansion for fermions (116) can be used to find the contribution to the effective action at zero temperature (vacuum energy). It rests on the following observation where E n is the energy of a mode with number n. It should be noted that the Fermi-Dirac distribution in this formula plays the role of a regularizing factor, and β 0 is a parameter of regularization. The derived formula gives the unrenormalized energy of vacuum fluctuations for one bosonic degree of freedom.
To get the vacuum energy, one can use the high-temperature fermionic expansion: Then for the vacuum energy we have Introduce the following counter-terms Consequently, the renormalized vacuum energy reads as which is twice as big as (antiparticle contribution included) the renormalized vacuum energy for a massive Let us sum the results up. We obtained the explicit expressions for the high-temperature expansion of the one-loop contribution to the thermodynamic potential of charged scalar particles in a constant crossed electromagnetic field. The vacuum energy was also calculated. The contributions of particles and antiparticles were separately investigated. It was shown explicitly that the high-temperature expansion of the thermodynamic potential and the vacuum energy do not contain contributions that are exponentially suppressed with respect to the external field or coupling constant except, possibly, the term at 1/β. The fact that the non-perturbative corrections to the vacuum energy do not depend on the external field for the given configuration of fields and plates can be anticipated from a simple analysis of the gauge invariant Lorentz-invariants. For the crossed electromagnetic field and given boundary conditions, the stress tensor has the form F µν := E 2 h [µ n ν] , where h µ = (1, 1, 0, 0), n µ = (0, 0, 0, 1) ⇒ hn = 0, n 2 = −1, h 2 = 0.
Here n µ is the normal to the hypersurface z = ±L/2. Due to antisymmetry of F µν , the following scalars vanish: η µν F µν = n µ n ν F µν = 0. All the higher powers of F µν are also zero, F k µν = 0, apart from F 2 µν = E 2 h µ h ν . Nevertheless, η µν F 2 µν = n µ n ν F 2 µν = 0. The same considerations apply to the dual tensorF µν = 1 2 ε µναβ F αβ , The cross products also vanish, (F F ) µν = 0. Notice that the presence of non-perturbative contributions in σ −1 ν that depend on the external field is not ruled out. These contributions may depend on the 4-vector specifying the reference frame where the thermodynamic system is at rest.
The system considered in the present paper can be used for an approximate description of thermodynamical properties of an ultrarelativistic fluid of charged bosons under the assumption that it is in a local thermodynamic equilibrium. In the reference frame comoving with the fluid element, the external electromagnetic field is crossed with good accuracy. Therefore, if the acceleration of this fluid element is small such that the effect of inertial forces on the thermodynamic properties of the system is negligible, then the thermodynamic properties of such a fluid element can be described by formulas obtained in the present paper.

A Dyson series
A heat kernel (the exponent on the right-hand side of (2)) can be considered as an evolution operator specified by the Hamiltonian −H(ω) (7) taken at imaginary time. The dependence of the evolution operator on the constant p 2 x + p 2 y + m 2 − ω 2 is trivial and will be taken into account only in the final answer. We use the perturbation theory (see, e.g., [2]) to find the matrix elements of the evolution operator taking H 0 =p 2 z as an unperturbed Hamiltonian and −2(ω + p x )Eẑ as a perturbation.
The solution to the Sturm-Liouville problem for the Hamiltonian H 0 with zero boundary conditions (7) is given by the following functions It is easy to see that in the basis formed by these functions the operator of coordinate and the unperturbed Hamiltonian have the following matrix elements: It should be noted that despite the fact that H 0 is quadratic in momentum [x, [x, [x, H 0 ]]] = 0. This is the reason why the naive expression for the heat kernel in external homogeneous electromagnetic fields [5,16] is not applicable in our case. Time evolution operator takes an especially simple form in the interaction picture (Dyson series) Consider the correction of the type E 0 : Taking into account the constant term discarded in −H(ω) and changing s = iτ , we get the zeroth approximation for the heat kernel where we set Re τ < 0 for the sake of convergence.
There is a trace of the heat kernel in the expression for the high-temperature expansion Sp G(ω, τ ) = S dp x dp y (2π) 2 G(ω, τ ; p x , p y ), where S := L x L y . The integrals over impulses are Gaussian and can be computed easily It's easy to find an asymptotic expansion of the remaining sum at small values of τ : In order to calculate the expansion of ζ + (ν, ω) at ω → ∞, it is convenient to use the following This gives the following contribution into the coefficients ζ k : There is no correction of the form E 1 owing to z nn = 0. Consider the correction of the form E 2 (denote it by index (2)): ∞ n=1 n|e −isĤ |n (2) ≈ −4(ω + p x ) 2 E 2 ∞ n,p=1 is L 2 π 2 z 2 np n 2 − p 2 e −is π 2 n 2 L 2 .

B Some series
Let us consider the series (the Epstein-Hurwitz zeta function, see [18,19]) that appears in the expression for σ 0 ν (µ) in (81). This series is a holomorphic function in the indicated domain of the complex α plane and can be continued by analyticity to the whole complex plane except a countable number of poles. Applying the Poisson summation formula to the series of the form (140), but with infinite limits, it is not difficult to derive where the contour C runs upward parallel to the imaginary axis such that Re(s + 2α) < −1 and Re s < −1.
The last two expressions in (142) can be used to construct the analytic continuation with respect to α to the region Re α ≤ −1. The last formula in (142) follows from the Mellin representation (see, e.g., [20]): x −1 ds s where the contour C goes upward from below parallel to the imaginary axis and to the left from the point s = −1.
The series, where x := π 2 /(4c 2 ). The latter series is reduced to the sum of zeta functions. The former series can be completed to the series with infinite limits, the term with n = 0 being understood as a limit. Then, using the Poisson formula, we deduce F 4 α (c) = ζ(4) + αζ (2) π 2 4c 2 − α(α − 1) 64 and the contour C runs upward parallel to the imaginary axis such that Re(s + 2α) < 3 and Re s < −1. The contributions (142), (147) are exponentially suppressed for c → +∞,