Algebraic Approach to Bose–Einstein Condensation in Relativistic Quantum Field Theory: Spontaneous Symmetry Breaking and the Goldstone Theorem

Abstract

We construct states describing Bose–Einstein condensates at finite temperature for a relativistic massive complex scalar field with $|\varphi |^4$-interaction. We start with the linearized theory over a classical condensate and construct interacting fields by perturbation theory. Using the concept of thermal masses, equilibrium states at finite temperature can be constructed by the methods developed in Fredenhagen and Lindner (Commun Math Phys 332:895, 2014) and Drago et al. (Ann Henri Poincaré 18:807, 2017). Here, the principle of perturbative agreement plays a crucial role. The apparent conflict with Goldstone’s theorem is resolved by the fact that the linearized theory breaks the U(1) symmetry; hence, the theorem applies only to the full series but not to the truncations at finite order which therefore can be free of infrared divergences.

1 Introduction

In this paper, we shall analyze the perturbative construction of a Bose–Einstein condensate for a relativistic charged scalar field theory at finite temperature.

The first experimental realization of Bose–Einstein condensation (BEC) in dilute vapors of alkali atoms has been obtained few years ago [3, 9, 22]. These works have pushed a lot the theoretical and experimental investigations of this phenomenon.

Usually, Bose–Einstein condensation is discussed in the realm of nonrelativistic quantum theories. (See, e.g., [58] and references therein.) Bose–Einstein condensation in the noninteracting case is the phenomenon that below a certain critical temperature the ground state becomes macroscopically populated. A similar phenomenon is seen in systems with interaction where an analogous interpretation in terms of an effective single-particle Hamiltonian can be made. The concept of BEC becomes mathematically precise by requiring that for an equilibrium state with average particle number N the largest eigenvalue of its one-particle reduced density matrix is at least of order N in the limit of large N. The eigenfunction to the highest eigenvalue of the reduced density matrix is then called the wave function of the condensate.

A basic fact about sufficiently dilute Bose gases is that the energy density is proportional to the square of the particle density, the proportionality factor being essentially the scattering length of the two-particle interaction potential, cf. Chapter 2 in [52]. By scaling the interaction potential so that the scattering length is proportional to the inverse of the particle number, one arrives in the large N limit at the Gross–Pitaevskii equation [35, 57] for the wave function of the condensate for Bosons that are confined in a trap, for instance, a finite box [49, 50, 52]. This particular limit (“GP-limit”) is different from the thermodynamic limit where the density as well as the interaction potential is not scaled with N.

Instead of considering a quantum system in a box and the limit where the box tends to cover the full space, we deal in this paper directly with the full spacetime. We shall work in a relativistic quantum field theory setting. This presents some advantages; first of all the construction of the algebra of interacting observables in this case is free from infrared divergences¹ thanks to the causal property of the theory, see, e.g., [12, 19].² Furthermore, there are physical systems where the relativistic nature of fundamental physics is manifest together with the phenomenon of condensation. Here we have in mind the possible existence of Boson stars at cosmological level [8] and condensation phenomena in high energy physics as, e.g., the quark gluon plasma [6, 62].

More precisely, we shall analyze various possible equilibrium states at finite temperature for a complex scalar quantum field with mass m and chemical potential $\mu$ both for free and self-interacting theories. In the free theory, if the mass m is larger than $|\mu |$, there is a single equilibrium state for a given temperature with n-point functions described by tempered distributions. If the mass m is smaller than $|\mu |$, there are no such states while if m equals $|\mu |$, there are various equilibrium states. These various states correspond to different phases of the system. Furthermore, the pure phases differ by macroscopic contributions; they have different one-point functions exhibiting spontaneous breakdown of the global U(1) symmetry. The charge density gets a finite contribution from the one-point functions. In the nonrelativistic limit, the charged scalar field tends to the nonrelativistic scalar field and the states tend to the known equilibrium states of a nonrelativistic system of spinless noninteracting bosons in the thermodynamic limit. The charge density tends to the particle density, and the nontrivial one-point function shows up in the long-distance behavior of the 2-point function which coincides with the one-particle-reduced density matrix in the thermodynamic limit. A discussion about the equivalence of BEC with spontaneous breaking of gauge symmetry in the nonrelativistic setting can be found, e.g., in [29, 51, 68].

We then consider thermal equilibrium states in the case of a $\varphi ^4$ self-interaction with positive coupling and chemical potential $\mu$. The traditional construction of nonzero temperature equilibrium states (KMS states) for perturbatively defined interacting quantum field theories suffers, even in the massive case, from spurious infrared divergences at higher-loop order (see [2, 64]). It was recently shown how these divergences can be circumvented [30, 53]. This construction amounts to an adaptation of formulas derived in rigorous statistical mechanics by Araki [4] to the framework of perturbative QFT. For $|\mu |<m$ the method works, and one obtains states which are invariant under the U(1)-symmetry.

For $|\mu |\ge m$, one expects spontaneous breakdown of symmetry. We thus expand the theory around a nontrivial solution $\phi$ of the classical field equation. In the limit of vanishing coupling constant $\lambda$ keeping the ratio $(\mu ^2-m^2)/\lambda$ finite, the classical background tends to the condensate of the free theory. If one instead keeps $\mu$ and m fixed, the classical background dominates over the quantum fluctuations; actually $|\phi |^2$ diverges as $\lambda ^{-1}$. In particular, if one scales the charge density by multiplying with $\lambda$, one gets the charge density of the classical solution rescaled by $\sqrt{\lambda }$ which is a solution of the field equation with $\lambda =1$. This limit can thus be seen as an analogous of the GP limit for relativistic systems.

Due to the Goldstone theorem, one has to deal with massless modes and therefore with a slow decay of correlation functions. In the case of a massless scalar field, this problem could be circumvented by taking into account that the interaction produces at finite temperature a thermal mass. If this term is included in the free theory, the correlations of the unperturbed state decay sufficiently fast [24]. In the case of BEC, this is in conflict with the existence of a Goldstone mode induced by the spontaneous breakdown of the U(1) symmetry of the model.

We solve this problem in the following way: We linearize the theory around the classical solution. The linearized theory breaks the U(1) symmetry and shows nonvanishing thermal masses; hence, the perturbative construction works as in the massless model treated in [24]. The U(1) symmetry is recovered for the full theory which then has a massless mode in agreement with the Goldstone theorem.

As mentioned before, we shall work in the relativistic quantum field theory setting called perturbative algebraic quantum field theory (pAQFT) [13, 14, 16, 32, 40, 41], see also the recent books on the subject [25, 60].

The paper is organized as follows: In the next section, we shall briefly recall the framework of pAQFT, the key steps in the construction of the KMS state performed in [30] and few facts about the principle of perturbative agreement [42] by which we can move the thermal mass into the free theory [24]. The third section contains the perturbative analysis of the massive complex scalar field with a $|\varphi |^4$ interaction, expanded around a solution with a nonvanishing condensate. We shall then discuss the construction of the interacting state at finite temperature over the condensate, and we analyze its adiabatic limit. The forth section contains the discussion of the formation of the condensate in connection to the spontaneous symmetry breaking. We shall actually see that the symmetry is effectively broken in the background theory while it is recovered in the exact theory where it is thus spontaneously broken.

2 Equilibrium States for Interacting Quantum Field Theory

In this section, we briefly review the formalism of perturbative algebraic quantum field theory. This framework combines renormalized perturbation theory with concepts of algebraic quantum field theory [36, 38]. The basic step is an assignment of the algebra generated by the observables of the theory to each region of the spacetime. Physical information, like locality, is then stored in the relations between algebras labeled by different regions of the spacetime. In the case of an interacting theory treated with perturbation theory, the elements of the algebras associated with each region are given as formal power series in the coupling constant with values in suitable $*$-algebras. States are constructed in a second step as linear functionals on the algebra of observables which via the GNS construction then provide representations of the elements of the algebra by operators on a state space. Renormalization in this framework is automatically independent of the state; moreover, infrared problems do not occur in the construction of the algebra. They may become visible in the construction of states where they indicate physical properties of the system.

2.1 Perturbative Construction of the Interacting Quantum Field Theory and the Adiabatic Limit

We shall here briefly recall the basic elements of the perturbative construction of the $\phi ^4$ scalar field theory propagating on a four dimensional Minkowski spacetime $({\mathbb {M}},\eta )$ where the metric $\eta$ has the signature $(-,+,+,+)$. The Lagrangian is

$$\begin{aligned} {\mathcal {L}} = -\frac{1}{2} \partial _\mu \phi \partial ^\mu \phi - \frac{m^2}{2}\phi ^2 - \frac{\lambda }{4} \phi ^4 \end{aligned}$$

(1)

where m is a positive mass and $\lambda$ the coupling to the nonlinear perturbation.

In a first step, we construct the algebra corresponding to the free theory ($\lambda =0$). We label the elements O of the algebra by functionals on the classical configuration space which in our case is the space of smooth functions $\phi$ on Minkowski space,

$$\begin{aligned} O[\phi ]= \sum _{n=0}^N \int f_n(x_1,\ldots , x_n) \phi (x_1)\ldots \phi (x_n)\mathrm{d}^{4n} x \end{aligned}$$

where $f_n$ is a compactly supported distribution on ${\mathbb {M}}^n$ which is symmetric under permutations of the arguments and where we used the measure induced by the Minkowski metric. The spacetime support $\mathrm {supp}\,O$ of O is the smallest closed subset G of Minkowski space such that $\text {supp}f_n\subset G^n$ for all n. O is called regular if all $f_n$ are smooth, and it is local if all $f_n$ are of the form

$$\begin{aligned} f_n(x_1,\ldots ,x_n)=P\prod _{i=2}^n\delta (x_1-x_i) \end{aligned}$$

with a partial differential operator P with smooth compactly supported coefficients.

The product of regular observables is defined in terms of the commutator function $\Delta$, which is the retarded minus advanced fundamental solution of $(\square -m^2)\phi =f$. In the sense of generating functionals, it is given by

$$\begin{aligned} e^{i\phi (f)}\star e^{i\phi (g)}=e^{-\frac{i}{2}\langle f,\Delta g\rangle }e^{i\phi (f+g)} \end{aligned}$$

where $\phi (f)=\int \phi (x)f(x)\mathrm{d}^4x$ for a compactly supported test function f. Due to singularities of $\Delta$, the product cannot be extended to nonlinear local functionals as, e.g., $\int f(x)\phi (x)^n\mathrm{d}^4x$ with $n>1$ and a test function f. This well-known problem is, as usual, circumvented by replacing these functionals by so-called normal-ordered functionals. In our framework, this means the following.

Consider the Klein–Gordon equation with spacetime-dependent mass m(x) on a globally hyperbolic spacetime. Let H be a symmetric bisolution of the form

$$\begin{aligned} H +\frac{i}{2}\Delta = \lim _{\epsilon \rightarrow 0^+}\frac{U}{\sigma _\epsilon } + V\log \left( \frac{\sigma _\epsilon }{\xi ^2}\right) +W \end{aligned}$$

(2)

where U, V and W are smooth symmetric functions of 2 spacetime points x and y. U and V depend only on the geometry and on the mass near to the geodesic connecting the arguments, and $\sigma _{\epsilon }(x,y) = \sigma (x,y) + i \epsilon (t(x) -t(y))$. $\sigma$ is the square of the geodesic distance between x and y, equipped with the appropriate sign for spacelike and timelike separation, respectively, t is a time function and $\xi$ a lengthscale. Such a bisolution is called a Hadamard function [47]. According to Radzikowski [59], a Hadamard function H can be characterized as a symmetric bisolution with the property that the wave front set of $H +\frac{i}{2}\Delta$ satisfies a positivity condition (microlocal spectrum condition [15]). Examples of Hadamard functions are the symmetric parts of the 2-point functions of vacuum and KMS states.

Given a Hadamard function H, normal ordering is a linear map defined by

$$\begin{aligned} e^{i\phi (f)}\mapsto :e^{i\phi (f)}:_H=e^{i\phi (f)}e^{\frac{1}{2}\langle f,Hf\rangle }. \end{aligned}$$

The $\star$-product then is interwined with the so-called Wick product $\star _H$,

$$\begin{aligned} :e^{i\phi (f)}:_H\star :e^{i\phi (g)}:_H=:e^{i\phi (f)}\star _H e^{i\phi (g)}:_H\equiv :e^{i\phi (f+g)}:_H e^{-\langle f,(H+\frac{i}{2}\Delta )g\rangle }. \end{aligned}$$

Due to the smaller wave front set of $H+\frac{i}{2}\Delta$ compared to $\Delta$, the Wick product can be extended to a larger class of functionals O. The product among these objects is well defined as long as $\text {WF}(f_n) \cap ({\overline{V}}_+^n\cup {\overline{V}}_-^n) = \emptyset$ where $\text {WF}(f_n)$ is the wave front set [43] of the distribution $f_n \in {\mathcal {E}}'({\mathbb {M}}^n)$ and where ${\overline{V}}_\pm ^n$ denotes the closure of the forward/backward light cones in $T^*{{\mathbb {M}}^n}$, namely ${\overline{V}}_\pm ^n \doteq \{(x_1,\ldots , x_n;k_1,\ldots , k_n)\in T^*M^n | \langle k_i, \eta ^{-1}k_i\rangle \le 0, k_i^0 \ge 0 \}$. These functionals are called microcausal, and we denote their set by ${\mathcal {F}}_{\mu c}$. It contains in particular the local functionals, denoted by ${\mathcal {F}}_{\text {loc}}$, and their pointwise (classical) products, the multilocal functionals. We refer to [13, 60] for further details on the definition of these sets. We can now extend the algebra of observables by normal-ordered microcausal functionals and define their product by

$$\begin{aligned} :O_1:_H\star :O_2:_H=:O_1\star _H O_2:_H\!\!. \end{aligned}$$

Note that the enlarged algebra ${\mathcal {A}}$, obtained as $({\mathcal {F}}_{\mu c},\star _H)$, does not depend on the choice of the Hadamard function. Only the labeling by the functional depends on H. Note, furthermore, that the normal-ordered functionals are, in general, no longer functionals on the configuration space due to the singularities of the Hadamard function.

The standard normal ordering in Fock space with respect to annihilation and creation operators is obtained if one chooses the symmetrized vacuum 2-point function $\Delta _1$ as the Hadamard function. It has the nice feature that the vacuum expectation value $\omega _0(:O:_{\Delta _1})$ of a normal-ordered functional $:O:_{\Delta _1}$ coincides with the evaluation of the functional O at $\phi =0$,

$$\begin{aligned} \omega _0(:O:_{\Delta _1})=O(\phi =0). \end{aligned}$$

A corresponding formula holds for any quasi-free (Gaussian) state $\omega$, if one uses its symmetrized 2-point function as the Hadamard function, in particular for KMS states.

This choice of normal ordering, however, is problematic when one wants to identify observables in different states or under the process of renormalization. Then, a preferred Hadamard function with $W=0$ is better behaved, as first discussed by Hollands and Wald [40]. In the case of a generic curved background that preferred function is in general not well defined, but in the case of Minkowski space, it is explicitly constructed in appendix A of [13]. It is unique up to the choice of a length scale.

The field equation is not yet implemented into the algebra ${\mathcal {A}}$. The algebra contains instead an ideal generated by normal-ordered functionals which vanish on solutions. The so-called on shell algebra is obtained by taking the quotient with respect to this ideal. This quotient is faithfully represented on Fock space and coincides with the standard algebra of the free field.

The off-shell algebra ${\mathcal {A}}$, however, is better behaved under the time-ordered product which is used for the incorporation of interaction.

Interacting fields can be constructed by means of causal perturbation theory, a method of renormalization elaborated by Epstein and Glaser [28] on the basis of ideas of Stueckelberg [66, 67] and Bogoliubov [7]. It was further developed by Scharf and collaborators (see, e.g., [63]). It is also the basis for a treatment of interactions on curved spacetimes [14, 40] where other versions of renormalization do not work. Its main idea is the construction of time-ordered products of interaction Lagrangians. In the work of Epstein and Glaser, these products are operator-valued distributions on Fock space. One uses the fact that the time-ordered product for noncoinciding points agrees with the operator product in the appropriate order. Renormalization then consists in extending these distributions to coinciding points. By induction with respect to the number of factors, one can show that this extension is always possible and unique up to the addition of a further interaction Lagrangian in each order. This corresponds precisely to the freedom in the choice of renormalization conditions known from other versions of renormalization.

In pAQFT, one uses a version of causal perturbation theory which is independent of the choice of a state space. There the time ordering operator is a linear map T from the algebra of multilocal functionals (with respect to the pointwise (classical) product)³ to the algebra of ${\mathcal {A}}$. On regular functionals, it is determined by

$$\begin{aligned} Te^{i\phi (f)}=e^{i\phi (f)}e^{-\frac{i}{2}\langle f,\Delta _D f\rangle }=:e^{i\phi (f)}e^{-\frac{i}{2}\langle f,\Delta _F^H f\rangle }:_H \end{aligned}$$

with the Dirac propagator $\Delta _D$ (the mean of retarded and advanced propagator) and the Feynman propagator $\Delta _F^H=\Delta _D-iH$ associated with the Hadamard distribution H. It satisfies the causal factorization condition

$$\begin{aligned} T(FG)=TF\star TG \end{aligned}$$

if the support of F does not intersect the past of the support of G. To extend T to multilocal functionals, one fixes a Hadamard function H and characterizes the extensions by the initial conditions

$$\begin{aligned} T1=1\text { and }TF=:F:_H \end{aligned}$$

for local functionals F vanishing at $\phi =0$. The causal factorization condition then fixes the map T on n-local functionals by its values on k-local functionals for $k<n$ up to a local functional.

In order to reduce the ambiguity in this extension, we choose a Hadamard distribution with $W=0$ as discussed before, so only the dependence on the scale $\xi$ remains for local functionals.

Given an interaction Lagrangian ${\mathcal {L}}_I$ where ${\mathcal {L}}_I(\varphi )$ is a translation invariant section of the jet bundle constructed over the field configuration $\varphi$, we consider the local functional

$$\begin{aligned} V[\varphi ] \doteq \lambda \int g {\mathcal {L}}_I(\varphi ) \mathrm{d}^4x \end{aligned}$$

(3)

with the test function g.

The formal S-matrix of the interaction Lagrangian V can now be constructed as time-ordered exponential

$$\begin{aligned} S(V) \doteq Te^{iV} \end{aligned}$$

in the sense of formal power series; hence, S(V) is an element of ${\mathcal {A}}[[\lambda ]]$, the set of formal power series in the coupling constant $\lambda$ with coefficients in ${\mathcal {A}}$. The $\star$-product in ${\mathcal {A}}$ extends directly to ${\mathcal {A}}[[\lambda ]]$. Relative S-matrices are then defined as

$$\begin{aligned} S_V(F) \doteq S(V)^{-1} \star S(V+F), \quad F\in {\mathcal {F}}_{\text {loc}} \end{aligned}$$

where the inverse is understood in terms of the $\star$ -product, and the interacting fields are given in terms of the Bogoliubov map (also called Møller operator) which extracts the contributions of $S_V(\mu F)$ linear in $\mu$

$$\begin{aligned} R_V(F) \doteq -i\left. \frac{\mathrm{d}}{\mathrm{d}\mu } S_V(\mu F) \right| _{\mu =0} = S(V)^{-1} \star T(e^{iV}F). \end{aligned}$$

Interacting observables can now be represented as elements of the algebra generated by the relative S matrices $S_V(F)$

$$\begin{aligned} {\mathcal {A}}_I({\mathcal {O}}) = \left[ \left\{ S_V(F) | F\in {\mathcal {F}}_{\text {loc}}({\mathcal {O}}) \right\} \right] , \end{aligned}$$

where the square brackets denote the set of linear combinations of products of the elements inside the brackets. We observe that the association of spacetime regions ${\mathcal {O}}$ to algebras ${\mathcal {A}}_I({\mathcal {O}})$ forms a net of subalgebras of ${\mathcal {A}}[[\lambda ]]$ in the sense of the Haag–Kastler axioms of algebraic quantum field theory. In the following, we omit the symbol $\star$ for the product within ${\mathcal {A}}[[\lambda ]]$ and replace it by juxtaposition.

The last step in constructing the interacting theory is the removal of the cutoff g from the interaction Lagrangian V in (3). This can be done taking the adiabatic limit $g\rightarrow 1$. At algebraic level, as soon as the interacting observables are supported on a compact region ${\mathcal {O}}$, this limit can be taken over larger and larger regions where the cutoffs are equal to 1; further details can be found in [14]. Here we can make this construction more explicit making use of the time slice axiom and the causal properties of the S-matrix. Actually, both the S-matrix and the relative S-matrix satisfy the following causal factorization property valid for $A,B,C\in {\mathcal {F}}_{\text {loc}}$

$$\begin{aligned} S(A+B+C) = S(A+B) S(B)^{-1} S(B+C), \quad A > rsim C \end{aligned}$$

(4)

where $A > rsim C$ means that A is later than C in the sense that $\text {supp}(A)\cap J_{-}(\text {supp}(C))= \emptyset$, where $J_\pm ({\mathcal {O}})$ denotes the causal future/past of ${\mathcal {O}}$. This causal factorization property implies that

$$\begin{aligned}&S_{V+W}(F) = S_V(F) , \quad W > rsim F \end{aligned}$$

(5)

$$\begin{aligned}&S_{V+W}(F) = S_V(W)^{-1} S_V(F) S_V(W) , \quad F > rsim W. \end{aligned}$$

(6)

If $g,g'$ coincide on $J_+({\mathcal {O}})\cap J_-({\mathcal {O}})$, the corresponding local functionals $V,V'$ differ by

$$\begin{aligned} V'-V=W_++W_- \end{aligned}$$

with $\text {supp}W_+\cap J_-({\mathcal {O}})=\emptyset$ and $\text {supp}W_-\cap J_+({\mathcal {O}})=\emptyset$; hence,

$$\begin{aligned} S_V(F)\mapsto S_{V'}(F)=S_V(W_-)S_V(F)S_V(W_-)^{-1}\text { for }F\in {\mathcal {F}}_{\text {loc}}({\mathcal {O}} \end{aligned}$$

extends to an isomorphism ${\mathcal {A}}^{g}_I({\mathcal {O}})\rightarrow {\mathcal {A}}^{g'}_I({\mathcal {O}})$. The limit $g\rightarrow 1$ can now be taken at the algebraic level.

Consider now a Cauchy surface $\Sigma = t^{-1}(0)$ where t is the time coordinate of a standard Minkowski coordinate system which is fixed once and forever. An $\epsilon$ neighborhood of the Cauchy surface $\Sigma$ is

$$\begin{aligned} \Sigma _\epsilon \doteq \{ p\in M | t(p)\in (-\epsilon ,\epsilon )\}. \end{aligned}$$

Interacting fields satisfy the time slice axiom, see [20]; namely, for every $A\in {\mathcal {A}}_{I}({\mathcal {O}})$ there exists a $C\in {\mathcal {A}}_{I}(\Sigma _\epsilon \cap J({\mathcal {O}}))$ such that

$$\begin{aligned} A=C+:W: \end{aligned}$$

where $W\in {\mathcal {F}}_{\mu c}$ vanishes on solutions $\phi$. Here $J({\mathcal {O}})=J_+({\mathcal {O}}) \cup J_-({\mathcal {O}})$. Hence, the on-shell algebras ${\mathcal {B}}_I({\mathcal {O}})$, obtained by taking quotients with respect to the ideal generated by these elements, are subsets of ${\mathcal {B}}_I(\Sigma _\epsilon )$. Thus, to construct a state for the interacting algebra ${\mathcal {B}}_{I}(M)$ it suffices to construct it for ${\mathcal {B}}_{I}(\Sigma _\epsilon )$. We therefore choose a cutoff function g of the form

$$\begin{aligned} g(t,{\mathbf {x}}) = \chi (t) h({\mathbf {x}}) \end{aligned}$$

(7)

where now the time cutoff is realized by $\chi (t)$ which is a smooth function which is equal to 1 for $t>-\epsilon$ and 0 for $t\le -2\epsilon$. Furthermore, h is a space cutoff which is compactly supported on $\Sigma$. To obtain a state in the adiabatic limit, it is sufficient to consider the limit where h tends to 1 keeping fixed the time cutoff $\chi$.

2.2 Interacting KMS States and the Adiabatic Limit

Equilibrium states are characterized by the Kubo–Martin–Schwinger (KMS) condition, see [37]. This condition yields canonical Gibbs states when they are well defined, but remains meaningful also for infinitely extended systems where the Gibbs formula can no longer be used [36, 37].

We recall here the definition. A state $\omega$ over a C*-algebra ${\mathfrak {B}}$ satisfies the KMS condition with respect to the one parameter group of $*$-automorphisms $\tau _t$ at inverse temperature $\beta$ if for every $A,B\in {\mathfrak {B}}$ $\omega (A\tau _t B)$ is an analytic function for $\text {Im}(t)\in (0,\beta )$, continuous at the boundary and if

$$\begin{aligned} \omega (A\tau _{i\beta } B) = \omega (BA). \end{aligned}$$

A state which satisfies the KMS condition at inverse temperature $\beta$ is called $\beta$-KMS state. It is automatically invariant under $\tau _t$ and satisfies similar relations for n-point functions.

If one uses the concept of KMS states for general $*$-algebras, one has to enrich the definition by some of these properties. See, e.g., Definition 1 in [30] for an extended discussion.

Let $\tau _t$ denote the one parameter group of $*$-automorphisms of ${\mathcal {B}}$, the $*$-algebra of free fields, induced by the action of Minkowski spacetime translations

$$\begin{aligned} \tau _t(F)\doteq F_t,\quad F_t[\varphi ] \doteq F[\varphi _t], \quad \varphi _t(s,{\mathbf {x}})\doteq \varphi (s-t,{\mathbf {x}}). \end{aligned}$$

Consider now the following two-point function

$$\begin{aligned} \omega _2^\beta (f,g) \doteq \frac{1}{(2\pi )^3}\int _{{\mathbb {R}}^4} \overline{{\hat{f}}(p) }{\hat{g}}(p) \frac{1}{1-e^{-\beta p_0}}\delta (p^{2}+m^2) \mathrm{d}^4p \end{aligned}$$

(8)

The quasi-free state $\omega ^\beta$ constructed out of this two-point function for the free theory is a KMS state at inverse temperature $\beta$ with respect to time translations. This state is easily described by using normal ordering with respect to its symmetrization $H_{\beta }$ as

$$\begin{aligned} \omega ^\beta (:F:_{H_{\beta }})= F[0], \quad F\in {\mathcal {F}}_{\mu c}. \end{aligned}$$

(9)

The interacting time evolution $\tau _t^V$ in ${\mathcal {B}}_I({\mathcal {O}})$, the subalgebra of ${\mathcal {B}}$ generated by S(F) with F supported in ${\mathcal {O}}$ which is a representation of the algebra of interacting observables supported in ${\mathcal {O}}$ is such that

$$\begin{aligned} \tau _t^V(S_V(F)) \doteq S_V(F_t), \end{aligned}$$

whereas the free evolution is $\tau _t(S_V(F)) = S_{V_t}(F_t)$. To construct the KMS state for the interacting theory, we have to relate the free and interacting time evolution. The causal factorization property (4) implies that

$$\begin{aligned} \tau ^V_t({S_V(F)}) = S_V(V_t-V) \tau _t(S_V(F)) S_V(V_t-V)^{-1},\quad S_V(F)\in {\mathcal {B}}_I(\Sigma _{\epsilon }),\quad t\ge 0. \end{aligned}$$

The map $t\mapsto U(t)\doteq S_V(V_t-V)$ from positive real numbers to unitary elements of ${\mathcal {B}}_I$ defines a cocycle which intertwines the free and interacting time evolution. The cocycle relation and its infinitesimal generators are

$$\begin{aligned} U(t+s) = U(t)\tau _t U(s), \quad H_I \doteq -i\left. \frac{\mathrm{d}}{\mathrm{d}t} U(t)\right| _{t=0} , \end{aligned}$$

where, in the case of V as in (3) with g as in (7), it turns out that

$$\begin{aligned} H_I = \int h({\mathbf {x}}) {\mathcal {H}}_I({\mathbf {x}}) \mathrm{d}^3{\mathbf {x}}, \quad {\mathcal {H}}_I({\mathbf {x}}) \doteq \int {{\dot{\chi }}}(t) R_V(-{\mathcal {L}}_I(t,{\mathbf {x}})) \mathrm{d}t. \end{aligned}$$

Hence, $H_I$ and ${\mathcal {H}}_I$ play the role of the interacting Hamiltonian and the interacting Hamiltonian density. We stress that due to the smearing in time, ${\mathcal {H}}_I({\mathbf {x}})$ is a well-defined formal power series with coefficients contained within the algebra of the free theory.

For any $\beta$-KMS state $\omega ^\beta$ of the free theory, like the quasi-free state (9) constructed with the two-point function (8) , and any spatial cutoff described by the test function h on $\Sigma$, we obtain a $\beta$-KMS state of the theory with interaction $H_I(h)$ with respect to the evolution $\tau ^V_t$ observing that

$$\begin{aligned} t\mapsto \omega ^{\beta }(A U(t)) \end{aligned}$$

for every $A\in {\mathcal {B}}_I(\Sigma _\epsilon )$ can be analytically continued to $\text {Im}t\in [0,\beta ]$. Hence,

$$\begin{aligned} \omega ^{\beta , V}_h (A) \doteq \frac{\omega ^{\beta }(A U(i\beta ))}{\omega ^{\beta } (U(i\beta ))}, \quad A\in {\mathcal {B}}_I(\Sigma _\epsilon ) \end{aligned}$$

(10)

defines a $\beta$-KMS state with respect to $\tau ^V_t$, as proved in [30].

Furthermore, the expectation values in the state $\omega _h^{\beta ,V}$ can be computed by the following formula

$$\begin{aligned} \omega ^{\beta ,V}_h(A)=&\sum _{n}\int _{0\le u_1\le \cdots u_n\le \beta } d u_1\ldots d u_n\int _{{{\mathbb {R}}}^{3n}}\mathrm{d}^3{\mathbf {x}}_1\ldots \mathrm{d}^3{\mathbf {x}}_n h({\mathbf {x}}_1)\ldots h({\mathbf {x}}_n)\nonumber \\&\omega _T^{\beta }\left( A;\tau _{iu_{1}}({\mathcal {H}}_I({\mathbf {x}}_{1}));\ldots ;\tau _{iu_n}({\mathcal {H}}_I({\mathbf {x}}_n))\right) , \quad A\in {\mathcal {B}}_I(\Sigma _\epsilon ) \end{aligned}$$

(11)

Here $\omega _T^\beta$ denotes the truncated functional associated with $\omega ^\beta$.

As shown in [30], the limit $h\rightarrow 1$ can now be taken provided the truncated n-point functions decay sufficiently fast for large spatial separations. Furthermore, the obtained state does not depend on $\chi$ anymore. In this way, one obtains the correlation functions for an interacting field in thermal equilibrium in the case of a massive theory.

2.3 Principle of Perturbative Agreement: Massless Case

If the linearized theory is massless, the limit $h\rightarrow 1$ in (11) cannot be taken because the decay of the n-point function is too slow. However, in the case of a $\phi ^4$ theory, it is possible to use a similar construction [24]. The idea is to modify the splitting

$$\begin{aligned} {\mathcal {L}} = {\mathcal {L}}_0+{\mathcal {L}}_I, \quad {\mathcal {L}}_0\doteq -\frac{1}{2}\partial \phi \partial \phi , \quad {\mathcal {L}}_I \doteq -\lambda \frac{1}{4} \phi ^4 \end{aligned}$$

(12)

by adding an artificial mass M to the background theory and subtracting it in the interacting Lagrangian, namely

$$\begin{aligned} {\mathcal {L}} = {\mathcal {L}}_0'+{\mathcal {L}}_I', \quad {\mathcal {L}}_0'\doteq {\mathcal {L}}_0 - \frac{{M}^2}{2} \phi ^2, \quad {\mathcal {L}}_I' \doteq {\mathcal {L}}_I + \frac{{M}^2}{2} \phi ^2. \end{aligned}$$

(13)

Let H be the distinguished Hadamard function for $\square -M^2$ with a length scale $\xi$ and T a time ordering operator with $TF=:F:_H$. Let $H_{\beta }=H+W_{\beta }$ denote the symmetrized 2-point function of the $\beta$-KMS state as in (8) for the theory with the modified free Lagrangian ${\mathcal {L}}_0'$. Then,

$$\begin{aligned} T(\phi ^4)=:\phi ^4:_{H_{\beta }}+6W_\beta (0,0):\phi ^2:_{H_\beta }+6W_{\beta }(0,0)^2. \end{aligned}$$

(14)

We see that the interaction Hamiltonian density in the KMS-state contains a mass term with a positive coefficient as long as $M^2<M^2_\beta$ with the thermal mass

$$\begin{aligned} M^2_{\beta }=3\lambda W_{\beta }(0,0). \end{aligned}$$

Under this condition, the interaction Lagrangian remains convex and possesses a single stationary point at $\phi =0$. As discussed for example in [24], the thermal mass is

$$\begin{aligned} M_\beta ^2 = \lambda \left( c_M M^2+ \frac{1}{2\pi ^2}\int _{0}^\infty \frac{1}{e^{\beta \sqrt{p^2+M^2}}-1} \frac{p^2}{\sqrt{p^2+M^2}} dp \right) \end{aligned}$$

where $c_M = \frac{1}{8\pi ^2} \log (M\xi )$ is a renormalization constant and it depends on the length scales $\xi$ in (2). If $\xi M$ is equal to 1, then $M_\beta ^2$ vanishes in the limit $\beta \rightarrow \infty$ and $M_\beta ^2 = T^2/12+O(M^2)$. We finally observe that the theories constructed with the two different splittings are equivalent thanks to the principle of perturbative agreement, which has been shown to hold in [42], see also [24]. Further details about the validity of this principle are collected in “Appendix A.” We finally recall that if the interaction Lagrangian is quadratic in the field $Q=\int \delta m^2 \varphi ^2 \mathrm{d}^4x$ and if it corresponds to a perturbation of the mass m of the free theory to $\sqrt{m^2+\delta m^2}$, the equilibrium state constructed as in (10) is the KMS state at inverse temperature $\beta$ with perturbed mass. This last observation has been proved in Theorem 3 of [23], and it shows that perturbative agreement is compatible with the construction of equilibrium states discussed in [30].

3 Massive Complex Scalar Field

We discuss the equilibrium states at finite temperature $\beta$ with a nonzero chemical potential $\mu$. We first discuss the free theory ($\lambda =0$); afterward, we study the corresponding states for the interacting theory. We are interested in finding states which can be interpreted as exhibiting Bose–Einstein condensation. The traditional way of defining BEC with particle numbers and occupation numbers cannot be applied in relativistic quantum systems. Instead, we look at states with nonvanishing one-point functions, thus showing spontaneous breakdown of the internal U(1) symmetry of the theory. This symmetry is generated by a conserved current $J_{\mu }$. The charge density $J_0$,

$$\begin{aligned} J_0(x)\doteq -i:{\dot{\varphi }}^*\varphi -\varphi ^*{\dot{\varphi }}:_H, \end{aligned}$$

replaces the particle density of the nonrelativistic theory. The mean of the charge density then distinguishes between different phases.

In the free theory, at fixed inverse temperature $\beta$, there is a critical value for the mean charge density. Below this value, the pure phases correspond to unique gauge invariant states, with a chemical potential $\mu$ depending on the charge density and with $\mu ^2<m^2$. If the charge density is above this threshold, the chemical potential has to satisfy $\mu ^2=m^2$; the states corresponding to pure phases have nonvanishing one-point functions which are related by the gauge symmetry. (See [10] for the concept of chemical potential in an algebraic formulation.)

Due to the nonvanishing one-point function, the two-point function is not decaying at large separations. Similar nonvanishing long- distance correlations are the basis of the criterion for BEC in the nonrelativistic theory. There one says that the ground state has a macroscopic occupation if the one-particle density matrix smeared in both entries over a spatial box of dimension L grows at least as particle number N in the limit where $L\rightarrow \infty$ keeping $N/L^3$ finite, see [52] for a more extensive discussion.

3.1 Condensate in the Free Theory

Let us start discussing the condensate for a free massive complex scalar quantum field theory propagating in a Minkowski spacetime. Let us denote by $\varphi$ the associated field configuration. Its equilibrium states with inverse temperature $\beta >0$ and chemical potential $\mu$, $|\mu |<m$ are the states which satisfy the KMS condition with respect to the time evolution

$$\begin{aligned} \tau _{t,\mu }(\varphi (x))=\varphi (x+te_0)e^{it\mu } \end{aligned}$$

(15)

where $e_0$ denotes the unit vector in time direction. The theory possesses an internal U(1) symmetry which might be spontaneously broken in some of equilibrium states. Hence, an interesting observable to distinguish these states is the current density

$$\begin{aligned} J_0(f)\doteq \int J_0(x)f(x) \mathrm{d}^4x \doteq -i\int \left( :{\dot{\varphi }}^*\varphi -\varphi ^*{\dot{\varphi }}:_H \right) f \mathrm{d}^4x \end{aligned}$$

(16)

where this is seen as an element of ${\mathcal {A}}$. We observe that, in view of the symmetry of the Hadamard coefficients U, $V_i$ of $V=V_n\sigma ^n$ in (2), and $W_i$ of $W=W_n\sigma ^n = (\omega _0-H-i\Delta /2)$, $\omega _0$ being the two-point function of the vacuum state, we have that $:{\dot{\varphi }}^*\varphi -\varphi ^*{\dot{\varphi }}:_H = :{\dot{\varphi }}^*\varphi -\varphi ^*{\dot{\varphi }}:_{H+W}.$ Hence, $J_0(f)$ can be seen as the current normal-ordered with respect to vacuum state. The possible pure phases are thus characterized by the following proposition:

Proposition 3.1

For inverse temperature $\beta >0$ and chemical potential $|\mu |<m$, there exists an unique KMS state with respect to $\tau _{t,\mu }$ whose n-point functions are tempered distributions. This state, denoted by $\omega _{\beta ,\mu }$ is quasi-free and its two-point functions are

$$\begin{aligned} \omega _{\beta ,\mu }(\varphi ^*(x)\varphi (y))\doteq \frac{1}{(2\pi )^{3}}\int {\mathrm{d}^4p}\,\delta (p^2+m^2)\epsilon (p_0)e^{ip(x-y)}\frac{1}{1-e^{-\beta (p_0-\mu )}} \end{aligned}$$

(17)

and $\omega _{\beta ,\mu }(\varphi (x)\varphi (y)))=\omega _{\beta ,\mu }(\varphi ^*(x)\varphi ^*(y)))=0$. The charge density in this state is

$$\begin{aligned} \omega _{\beta ,\mu }(J_0(x))&=\int {d}^4p\,2|p_0|\delta (p^2+m^2)\left( \frac{1}{e^{\beta (|p_0|-\mu )}-1}-\frac{1}{e^{\beta (|p_0|+\mu )}-1}\right) . \end{aligned}$$

It holds that

$$\begin{aligned} |\omega _{\beta ,\mu }(J_0(x)) | \le \rho _{cr}(\beta )\doteq \omega _{\beta ,m}(J_0(x)) \end{aligned}$$

(18)

where $\rho _{cr}(\beta )$ is the critical charge density. For $\beta >0$ and $\mu = \pm m$, there exist various KMS states with respect to $\tau _{t,\mu }$. Let us denote by $\Omega _{\beta ,\pm m}$ the set of quasi-free KMS states. The pure phases are the extremal points in $\Omega _{\beta ,\pm m}$, and these states are

$$\begin{aligned} \omega _{\beta ,c}^\pm = \omega _{\beta ,\pm m}\circ \gamma _c^{\pm } \end{aligned}$$

where $\gamma _c^{\pm }$ is an automorphism which is generated by

$$\begin{aligned} \gamma ^{\pm }_c(\varphi (x))=\varphi (x)+e^{\pm ix^0m}c({\mathbf {x}}) \end{aligned}$$

(19)

where c is a harmonic function of the spatial variables ${\mathbf {x}}$, $\Delta c=0$, with the spatial Laplacian $\Delta$. In this case,

$$\begin{aligned} \omega _{\beta ,c}^{\pm }(J_0(x))&=\int {d}^4p\,2|p_0|\delta (p^2+m^2)\left( \frac{1}{e^{\beta (|p_0|\mp m)}-1}-\frac{1}{e^{\beta (|p_0|\pm m)}-1}\right) \pm 2m|c|^2. \end{aligned}$$

Furthermore,

$$\begin{aligned} |\omega _{\beta ,c}^{\pm }(J_0(x))| \ge \rho _{cr}(\beta ). \end{aligned}$$

Proof

First of all, we observe that the KMS states corresponding to pure phases are quasi-free states with at most a nontrivial one-point function. A proof of this fact can be found in [61]. Furthermore, the truncated two-point function is constrained by the KMS condition to be equal to (17). The one-point function $\omega (\varphi )$ is constrained by the equation of motion and by request of invariance under the action of $\tau _{t,\mu }$. In particular, invariance under the action $\tau _{t,\mu }$ implies that the function

$$\begin{aligned} t\mapsto \omega (\varphi (x))e^{-i\mu t} \end{aligned}$$

is constant in time. This function needs to be a solution of $-\Delta + m^2 -\mu ^2$; however, for $|\mu |<m$ these solutions cannot be tempered distributions. The inequalities involving the critical charge density $\rho _{cr}(\beta )$ are an immediate consequence of the form of the expectation value of $J_0$ in the analyzed states.

$\square$

The state $\omega _{\beta ,\mu }$ respects the U(1)-symmetry of the theory. For chemical potentials $\mu =\pm m$, there exist many equilibrium states at fixed temperature and the U(1) symmetry is spontaneously broken in the states $\omega _{\beta ,c}^{\pm }$ for $c\ne 0$. At zero temperature, all the states with chemical potential $|\mu |<m$ coincide with the vacuum. Hence, the vacuum expectation value of the charge density $J_0(x)$ vanishes in that limit. Thus, any nonvanishing charge density in the limit $T=\beta ^{-1}\rightarrow 0$ of vanishing temperature requires a condensate.

At finite $\beta$ and $\mu =m$, we have $\omega _{\beta ,c}^{+}(J_0(x)) = \rho _{cr}(\beta ) + 2m|c|^2$. A nonvanishing condensate can occur only if $|\omega _{\beta ,c}^{+}(J_0)|>\rho _{cr}(\beta )$. Since $\rho _{cr}(\beta )>0$ is monotonically decreasing in $\beta$, diverges for $\beta \rightarrow 0$ and tends to 0 for large $\beta$, at a fixed charge density $\omega _{\beta ,c}^{\pm }(J_0)$, there is a critical temperature $T_{cr}>0$ such that only for $T<T_{cr}$ a condensate can be formed.

In “Appendix B,” we shall compute the nonrelativistic limit of the states analyzed in this section showing that the charged scalar field tends to the nonrelativistic scalar field and the states tend to the known equilibrium states of a nonrelativistic system of spinless noninteracting bosons in the thermodynamic limit. Furthermore, the charge density converges to the particle density. Finally, we see that the nontrivial one-point function shows up in the long-distance behavior of the 2-point function which coincides with the one-particle-reduced density matrix in the thermodynamic limit.

3.2 Massive Complex Scalar Field with $\varphi ^4$ Interaction over the Condensate

In this section, we start discussing the perturbative construction of the $\varphi ^4$ interacting theory over a suitable classical solution of the equation of motion which represents the condensate in the Minkowski spacetime. The Lagrangian of the theory we are considering is thus

$$\begin{aligned} {\mathcal {L}} = -\frac{1}{2} \partial {\overline{\varphi }}\partial {\varphi } - \frac{1}{2}m^2 |\varphi |^2 - \frac{\lambda }{4} |\varphi |^4 \end{aligned}$$

where $\varphi$ is a complex scalar field. Following a similar procedure presented in section III of [1], we expand ${\mathcal {L}}$ around a real classical solution $\phi$ which represents the condensate. Hence,

$$\begin{aligned} \varphi = e^{-i\mu x^0}(\phi + \psi ) \end{aligned}$$

(20)

where $\mu$ is again the chemical potential, $x^0$ is a fixed Minkowski time and $\psi$ is a complex scalar field which describes the perturbations. Its real and imaginary parts are denoted by $\psi _1$ and $\psi _2$, and thus,

$$\begin{aligned} \psi = \psi _1+i\psi _2. \end{aligned}$$

The Lagrangian density can now be written as a sum of contributions homogenous in the number of fields $\psi$ as follows:

$$\begin{aligned} {\mathcal {L}} = {\mathcal {L}}_0+{\mathcal {L}}_2+{\mathcal {L}}_3+{\mathcal {L}}_4 \end{aligned}$$

where

$$\begin{aligned} {\mathcal {L}}_0&= \frac{1}{2} \left| (\partial _0-i\mu )\phi \right| ^2 - \frac{1}{2} |\nabla \phi |^2 - \frac{\lambda }{4}|\phi |^4 - \frac{1}{2}m^2|\phi |^2 \\ {\mathcal {L}}_2&= \frac{1}{2} \left| (\partial _0-i\mu )\psi \right| ^2 - \frac{1}{2} |\nabla \psi |^2 -\lambda \phi ^2|\psi _1|^2 -\frac{1}{2}(\lambda \phi ^2+m^2)|\psi |^2 \\ {\mathcal {L}}_3&= -\lambda \phi \psi _1 |\psi |^2 \\ {\mathcal {L}}_4&= -\frac{\lambda }{4}|\psi |^4. \end{aligned}$$

The term ${\mathcal {L}}_1$ vanishes because $\phi$ is chosen to be a stationary point for the classical action $\int {\mathcal {L}}_0 \mathrm{d}^4x$. In the following, we shall choose a nonvanishing $\phi$ to describe the condensate, we discuss the quantization of the linearized theory (${\mathcal {L}}_2$), and finally, we use perturbation theory over the linearized theory to take into account ${\mathcal {L}}_3+{\mathcal {L}}_4$.

Contrary to the case of the free theory, in the interacting theory the chemical potential is not restricted to the interval $[-m,m]$. A chemical potential outside of this interval induces a spontaneous breakdown of symmetry showing up in a nonvanishing one-point-function and, as a consequence, in long-range behavior of the two-point-function, similar to the nonrelativistic case. In contrast to the free case, states with different condensates are not in mutual thermal equilibrium, since their chemical potentials differ.

3.2.1 The Condensate in the Vacuum Theory

We look for the case of a translation invariant background $\phi$. Then, the kinetic term in ${\mathcal {L}}_0$ has no effect and $\phi$ is a stationary point for

$$\begin{aligned} I=\int U(|\phi ^2|) \mathrm{d}^4x \end{aligned}$$

where

$$\begin{aligned} U(|\phi ^2|) = -\frac{\lambda }{4}|\phi |^4-\frac{1}{2}(m^2-\mu ^2)|\phi |^2; \end{aligned}$$

hence, it holds

$$\begin{aligned} |\phi |^2 = \frac{\mu ^2-m^2}{\lambda } \end{aligned}$$

(21)

and only one real, positive and translational invariant background solution $\phi$ is thus available for $\mu ^2>m^2$. We notice that for fixed $\mu ^2>m^2$ the background value of the field $\phi$ is of order $1/\sqrt{\lambda }$. In this case, we observe that ${\mathcal {L}}_2$ does not depend on $\lambda$, ${\mathcal {L}}_3$ is of order $\sqrt{\lambda }$ while ${\mathcal {L}}_4$ is of order $\lambda$. In the next, we shall construct the interacting field theory with perturbation methods considering ${\mathcal {L}}^I=({\mathcal {L}}_3+{\mathcal {L}}_4)$ the interaction Lagrangian. Hence, the solution we shall obtain will be a formal power series in $\sqrt{\lambda }$.

In the next, we shall discuss the construction of the quantum theory over the background discussed so far.

We argue that there exists a limit in which all the correlation functions are dominated by the classical background $\phi$. Actually, in the limit $\lambda \rightarrow 0$ keeping $|\mu |^2-m^2$ finite, the classical background $\phi$ diverges as $\lambda ^{-1/2}$, furthermore, the linearized theory is not affected by changes of $\lambda$ while the S-matrix constructed with the interacting Lagrangian tends to 1.

Hence, we expect that, under this limit, the one-point function rescaled by $\sqrt{\lambda }$ tends to the background value ${\tilde{\phi }} e^{-i m x^0}$ where ${\tilde{\phi }} = \sqrt{\lambda } \phi$, and similarly, the rescaled charge density $\lambda J_0$ tends to the charge density of the background $2\mu |{\tilde{\phi }}|^2$. Both these quantities do not depend on $\lambda$. We finally observe that the rescaled background ${\tilde{\phi }}$ is a solution of the equation of motion descending from the rescaled classical Lagrangian density

$$\begin{aligned} \tilde{{\mathcal {L}}}_0 = \lambda {\mathcal {L}}_0 = \frac{1}{2} |(\partial _0-i\mu ){\tilde{\phi }}|^2 - \frac{1}{2} |\nabla {\tilde{\phi }}|^2 - \frac{1}{4}|{\tilde{\phi }}|^4 - \frac{1}{2}m^2|{\tilde{\phi }}|^2 \end{aligned}$$

which is also independent on $\lambda$.

This is in analogy to what happens in the nonrelativistic case; actually there under the Gross–Pitaevskii limit, the density of the ground state tends to the density of a suitable classical solution of the Gross–Pitaevskii equation [49, 50], see in particular Theorem 1.1 and Theorem 1.2 in [50].

We thus argue that the equation of motion corresponding to the rescaled zeroth-order Lagrangian $\tilde{{\mathcal {L}}}_0$ can be interpreted as an analogous of the Gross–Pitaevskii equation in the relativistic setting, and thus, the limit $\lambda \rightarrow 0$ taken with m and $\mu$ fixed can be understood as the analogous of the Gross–Pitaevskii limit discussed in introduction.

3.2.2 Linearized Theory

The first step to construct the quantization of $\varphi$ is the analysis of the linearized equations of motion for the fluctuations $(\psi _1,\psi _2)$ around $\phi$. They have the form

$$\begin{aligned} \begin{aligned} (\square -M_1^2)\psi _1-2\mu {\dot{\psi }}_2&=0 \\ (\square -M_2^2)\psi _2+2\mu {\dot{\psi }}_1&=0 \end{aligned} \end{aligned}$$

(22)

where

$$\begin{aligned} M_1^2=(m^2-\mu ^2) + 3\lambda \phi ^2 \quad \text {and}\quad M_2^2= (m^2-\mu ^2)+ \lambda \phi ^2. \end{aligned}$$

(23)

Notice that if (21) holds, $M_1^2=2(\mu ^2-m^2)$ and $M_2^2=0$. Hence, we assume $M_1>M_2\ge 0$. Let us introduce

$$\begin{aligned} M^2 =\frac{M_1^2+M_2^2}{2}, \quad \delta M^2 =\frac{M_1^2-M_2^2}{2}. \end{aligned}$$

We observe that Eq. (22) for ${\tilde{\psi }}=(\psi _1,\psi _2)$ can be written in a compact form $D {{\tilde{\psi }}} = 0$ where D is given in terms of the standard Pauli matrices

$$\begin{aligned} {{\varvec{\sigma }}}_1=\begin{pmatrix} 0&{}\quad 1\\ 1&{}\quad 0 \end{pmatrix}, \quad {{\varvec{\sigma }}}_2=\begin{pmatrix} 0&{}\quad -i\\ i&{}\quad 0 \end{pmatrix} ,\quad {{\varvec{\sigma }}}_3=\begin{pmatrix} 1&{}\quad 0\\ 0&{}\quad -1 \end{pmatrix}. \end{aligned}$$

as

$$\begin{aligned} \begin{aligned} D&= (\square -M^2) {\mathbb {I}} - \delta M^2 {{\varvec{\sigma }}}_3 - i 2\mu \partial _0 {{\varvec{\sigma }}}_2 \\ {\overline{D}}&= (\square -M^2) {\mathbb {I}} + \delta M^2 {{\varvec{\sigma }}}_3 + i 2\mu \partial _0 {{\varvec{\sigma }}}_2. \end{aligned} \end{aligned}$$

(24)

Notice that

$$\begin{aligned} D{\overline{D}}= \left( (-\square +M^2)^2-(\delta M^2)^2 + 4\mu ^2 \partial _0^2\right) {\mathbb {I}}. \end{aligned}$$

The retarded and advanced propagators of the theory can be obtained as

$$\begin{aligned} \Delta _R \doteq {\overline{D}} (D{\overline{D}})_R,\quad \Delta _A \doteq {\overline{D}} (D{\overline{D}})_A, \end{aligned}$$

where $(D{\overline{D}})_R$ and $(D{\overline{D}})_A$ are the retarded and advanced fundamental solutions of $D{\overline{D}}$. Let us thus study

$$\begin{aligned} \widehat{D{\overline{D}}}&= \left( (p^2+M^2)^2 - (\delta M^2)^2 -4\mu ^2p_0^2 \right) {\mathbb {I}} \\&= \left( p_0^4 -2p_0^2 ({\mathbf {p}}^2+M^2 +2\mu ^2) + ({\mathbf {p}}^2+M^2)^2 - (\delta M^2)^2\right) {\mathbb {I}} \end{aligned}$$

Hence, the four solutions of $p_0^4 -2p_0^2 ({\mathbf {p}}^2+M^2 +2\mu ^2) + ({\mathbf {p}}^2+M^2)^2 - (\delta M^2)^2=0$ are $\pm \omega _\pm$ where

$$\begin{aligned} \omega _\pm ^2&= w^2 +2\mu ^2 \pm \sqrt{(w^2 +2\mu ^2)^2 - w^4 +(\delta M^2)^2} \nonumber \\&= w^2 +2\mu ^2 \pm \sqrt{4\mu ^4 + 4\mu ^2w^2 +(\delta M^2)^2} \nonumber \\&= w^2 +2\mu ^2 \pm \sqrt{(w^2 +2\mu ^2)^2 - w_1^2w_2^2} \end{aligned}$$

(25)

where now $w^2 \doteq {\mathbf {p}}^2+M^2$ and $w_i^2\doteq {\mathbf {p}}^2+M^2_i$.

We notice that if $M_2=0$ we have that $w_2=0$ for $|{\mathbf {p}}|=0$, and thus,

$$\begin{aligned} \lim _{|{\mathbf {p}}| \rightarrow 0 }\omega _-^2 = 0; \end{aligned}$$

hence, a massless mode is present in this system as expected by the Goldstone theorem. However, if the linearized theory is not in a ground state, it could happen that the normal-ordered interaction Lagrangian with respect to the state, as in (14), contains quadratic terms that could contribute to the masses of the fluctuations. If we use the formula

$$\begin{aligned} \prod _i \frac{1}{x-x_i} = \sum _i \frac{1}{x-x_i} \prod _{j\ne i} \frac{1}{x_i-x_j} \end{aligned}$$

valid for pairwise different $x_1,\ldots , x_n$, a couple of times, we get

$$\begin{aligned} {{\widehat{\Delta }}}_R&= \frac{\widehat{{\overline{D}}}}{(\omega _+^2-\omega _-^2)} \left( \frac{1}{(p_0+i\epsilon )^2-\omega _+^2} -\frac{1}{(p_0+i\epsilon )^2-\omega _-^2}\right) \end{aligned}$$

where recalling (24)

$$\begin{aligned} \widehat{{\overline{D}}} = -(p^2+M^2) {\mathbb {I}} + \delta M^2 {{\varvec{\sigma }}}_3 + 2\mu p_0 {{\varvec{\sigma }}}_2. \end{aligned}$$

(26)

We can construct $\Delta _A$ just changing $i\epsilon \rightarrow - i\epsilon$, while the Feynman propagator $\Delta _F$ is obtained substituting $(p_0+i\epsilon )^2$ with $p_0^2+i\epsilon$ and multiplying by i. Finally, the commutator function is

$$\begin{aligned} {{\widehat{\Delta }}} = \frac{2\pi i \widehat{{\overline{D}}}}{\omega _+^2-\omega _-^2} \epsilon {(p_0)}\left( \delta (p_0^2-\omega _+^2)-\delta (p_0^2-\omega _-^2) \right) \!. \end{aligned}$$

With $\Delta$ at disposal, the quantum product can be given as in Sect. 2.1; in this way we obtain the $*$-algebra of field observables ${\mathcal {A}}_0$. The analog of the Hadamard singularity H (2) for this theory can be given. The form of some of the corresponding Hadamard coefficients is discussed in “Appendix C.” The extended $*$-algebra of field observables ${\mathcal {A}}$ containing Wick polynomials normal-ordered with respect to H is obtained as in 2.1.

3.3 KMS States for the Linearized Theory

In view of the decomposition of the field $\varphi$ given in (20), the action of $\tau _t$ on $\psi$ as time translation is equivalent to the action of $\tau _{t,\mu }$ on $\varphi$ as given in (15). Hence, having the causal propagator of the linearized theory at disposal, we can construct the two-point function of the quasi-free $\beta$-KMS state with respect to time translation $\tau _t$ of the $\psi$ fields as

$$\begin{aligned} {\widehat{\omega }}_{\beta ,\psi } = \frac{i{{\widehat{\Delta }}}}{1-e^{-\beta p_0}}. \end{aligned}$$

Introducing

$$\begin{aligned} {\mathcal {S}} \doteq \begin{pmatrix} \psi _1(x)\psi _1(y) &{}\quad \psi _1(x)\psi _2(y)\\ \psi _2(x)\psi _1(y)&{}\quad \psi _2(x)\psi _2(y) \end{pmatrix} \end{aligned}$$

we have that the two-point function of the quasi-free $\beta$-KMS state $\omega _{\beta ,\psi }$ is in position space

$$\begin{aligned} \omega _{\beta ,\psi }({\mathcal {S}})= & {} \frac{1}{(2\pi )^{3}} \int {d}^4p\; e^{ip(x-y)} \frac{ \epsilon {(p_0)}}{\omega _+^2-\omega _-^2} \nonumber \\&\left( \delta (p_0^2-\omega _+^2)-\delta (p_0^2-\omega _-^2) \right) \frac{(-\widehat{{\overline{D}}})}{1-e^{-\beta p_0}}. \end{aligned}$$

(27)

Recalling the form of $\omega _\pm$ in (25), we notice that if $M_1>M_2>0$

$$\begin{aligned}&\omega _\pm ^2=w^2 +2\mu ^2 \pm \sqrt{(w^2 +2\mu ^2)^2 - w_1^2w_2^2}>0, \\&\omega _+^2-\omega _-^2 =2\sqrt{4\mu ^4 + 4\mu ^2w^2 +(\delta M^2)^2} >0, \end{aligned}$$

this means that no infrared divergences are present in $\omega _{\beta ,\psi }$ if $M_2>0$. The two-point function of the ground state of the $\psi _i$ theory (keeping the condensate $\phi \ne 0$) can be obtained taking the limit $\beta \rightarrow \infty$ of (27). Hence, to study expectation values in the state $\omega _{\beta ,\psi }$ of observables normal-ordered with respect to the vacuum $\omega _{\infty ,\psi }$ we consider $W=\omega _{\beta ,\psi }-\omega _{\infty ,\psi }$ and we obtain

$$\begin{aligned} W({\mathcal {S}})&= \frac{1}{(2\pi )^{3}} \int {d}^4p\; e^{ip(x-y)} \frac{1}{\omega _+^2-\omega _-^2} \left( \delta (p_0^2-\omega _+^2)-\delta (p_0^2-\omega _-^2) \right) \frac{(-\widehat{{\overline{D}}})}{e^{\beta |p_0|}-1}. \end{aligned}$$

(28)

We observe that in the coinciding point limit, the off-diagonal expectation values are vanishing

$$\begin{aligned} W(\psi _1(x)\psi _2(x))=0, \quad W(\psi _2(x)\psi _1(x))=0, \end{aligned}$$

and introducing $2\delta \omega ^2\doteq \omega _+^2-\omega _-^2 = 2\sqrt{4\mu ^4 + 4\mu ^2 w^2 +(\delta M^2)^2}$ we have that $W(\psi _i^2)$; the coinciding point limits of the diagonal elements of $W({\mathcal {S}})$ are

$$\begin{aligned} W(\psi _1^2)&= \frac{1}{(2\pi )^{3}} \int {d}^3{\mathbf {p}}\; \left( \frac{\delta \omega ^2+2\mu ^2+\delta M^2 }{\delta \omega ^2} \frac{1}{2\omega _+} \frac{1}{e^{\beta \omega _+}-1} \right. \nonumber \\&\quad \left. + \frac{ \delta \omega ^2 -2\mu ^2-\delta M^2 }{\delta \omega ^2}\frac{1}{2\omega _-}\frac{1}{e^{\beta \omega _-}-1} \right) \end{aligned}$$

(29)

$$\begin{aligned} W(\psi _2^2)&= \frac{1}{(2\pi )^{3}} \int {d}^3{\mathbf {p}}\; \left( \frac{\delta \omega ^2+2\mu ^2-\delta M^2 }{\delta \omega ^2} \frac{1}{2\omega _+} \frac{1}{e^{\beta \omega _+}-1} \right. \nonumber \\&\quad \left. + \frac{ \delta \omega ^2 -2\mu ^2+\delta M^2 }{\delta \omega ^2}\frac{1}{2\omega _-}\frac{1}{e^{\beta \omega _-}-1} \right) \end{aligned}$$

(30)

Notice that the integrand in both $W(\psi _1^2)$ and $W(\psi _2^2)$ is positive.

To analyze some properties of the condensate in the linearized theory, we compare the expectation values of the current density in the state $\omega _{\beta ,\psi }$ with (16), namely the current density of the free theory analyzed in Sect. 3.1. To this end, we recall the decomposition (20) and we get that

$$\begin{aligned} J_0=-i (:\dot{{\overline{\varphi }}}{\varphi }-{{\overline{\varphi }}}{\dot{\varphi }}:_H) = {\tilde{j}} -i \phi (\dot{{\overline{\psi }}}-{\dot{\psi }})+2\mu :|\phi +\psi |^2 :_H \end{aligned}$$

(31)

where now

$$\begin{aligned} {\tilde{j}}=-i \left( :\dot{{\overline{\psi }}}\psi -{\overline{\psi }}{{\dot{\psi }}} :_H\right) = 2\left( :{{\dot{\psi }}}_1\psi _2-\psi _1{{\dot{\psi }}}_2:_H\right) \end{aligned}$$

and H is the distinguished Hadamard function constructed in (57). We furthermore observe that, up to some choice of the renormalization freedom, $:{\dot{\psi }}_1\psi _2-\psi _1{\dot{\psi }}_2:_H=:{\dot{\psi }}_1\psi _2-\psi _1{\dot{\psi }}_2:_{\omega _{\infty ,\psi }}$ and $:|\psi |^2:_{H}=:|\psi |^2:_{\omega _{\infty ,\psi }}=:\psi _1^2:_{\omega _{\infty ,\psi }}+:\psi _2^2:_{\omega _{\infty ,\psi }}$. Hence,

$$\begin{aligned} \omega _{\beta ,\psi }(:|\psi |^2:_H)&=\frac{1}{(2\pi )^{3}} \int {d}^3{\mathbf {p}} \nonumber \\&\quad \left( \left( 1+ \frac{2\mu ^2 }{\delta \omega ^2}\right) \frac{1}{\omega _+} \frac{1}{e^{\beta \omega _+}-1}+ \left( 1-\frac{2\mu ^2}{\delta \omega ^2}\right) \frac{1}{\omega _-}\frac{1}{e^{\beta \omega _-}-1} \right) \end{aligned}$$

(32)

and

$$\begin{aligned} \omega _{\beta ,\psi }({\tilde{j}})&= \frac{4\mu }{(2\pi )^{3}} \int {d}^3{\mathbf {p}}\; \frac{1}{\delta \omega ^2} \left( \frac{\omega _-}{e^{\beta \omega _-}-1} - \frac{\omega _+}{e^{\beta \omega _+}-1} \right) . \end{aligned}$$

(33)

Hence,

$$\begin{aligned} \omega _{\beta ,\psi }(J_0) = \omega _{\beta ,\psi }({\tilde{j}}) +2\mu \; \omega _{\beta ,\psi }(:|\psi |^2:_H) +2\mu |\phi |^2. \end{aligned}$$

Notice that $\omega _+^2>\omega _-^2$ and that $\delta \omega ^2 \ge 2 \mu ^2$; hence, the integrand in $\omega _{\beta ,\psi }(:|\psi |^2:_H)$ given in (32) is always positive and monotonically decreasing in $\beta$. Similarly, for positive $\mu$, the integrand in $\omega _{\beta ,\psi }({\tilde{j}})$ given in (33) is also always positive and monotonically decreasing in $\beta$. Finally, both expressions (32) and (33) are diverging for $\beta \rightarrow 0$ and vanishes for $\beta \rightarrow \infty$. Hence, similar to the discussion given in Sect. 3.1 we have that $2\mu |\phi ^2|$ plays the role of the condensate charge density.

Consider now the case where $M_i$ are given in (23) with $\phi$ chosen to satisfy (21). In this case, $\lambda |\phi |^2 = \mu ^2-m^2$, and thus, the linearized theory does not depend on $\lambda$ while the background field scales as $\lambda ^{-1}$. The charge density is thus dominated by the charge density of the background $2\mu |\phi |^2$ thus confirming that the limit $\lambda \rightarrow 0$ taken with fixed $\mu$ and m is the relativistic analogous of the Gross–Pitaevskii limit discussed in introduction and at the end of Sect. 3.2.1.

Following closely the discussion given at the end of Sect. 3.1, we also see that in this case the critical charge density equals $\rho _{cr}(\beta )$ given implicitly in (18). Finally, in the limit $\lambda \rightarrow 0$ taken keeping the ratio $(\mu ^2-m^2)/\lambda$ finite we have that the states $\omega _{\beta ,\psi }$ of the linearized theory discussed so far tend to $\omega ^\pm _{\beta ,c}$ with $c=\phi$.

3.3.1 Thermal Masses

Having analyzed the equilibrium state of the free theory on field observables ${\mathcal {A}}$, the next step in the construction of an equilibrium state for the interacting theory will be an application of the analysis given in [30] and summarized in Sect. 2.2, namely to use (10) starting with a quasi-free state whose two-point function is given in (27). However, we expect that the limit $h\rightarrow 1$ cannot be directly taken because, as discussed above, if (21) holds, the mass $M_2$ given in (23) vanishes; hence, for vanishing spatial momentum, $\omega _-^2$ is also vanishing. This implies that various propagators of the linearized theory diverge for $p\rightarrow 0$. Hence, in agreement with Goldstone theorem a massless mode is present in this case. This implies a slow decay in the connected n-point functions constructed with $\omega _{\beta , \mu }$ given in (27).

In order to cure this problem, we use a different splitting of the Lagrangian into the free and interacting part. Actually, we add a virtual mass $m_v^2$ to the linearized fields and we remove them in the interaction Lagrangian. More precisely, the Lagrangian of the free theory is now

$$\begin{aligned} {\mathcal {L}}_2' = {\mathcal {L}}_2 - \frac{m_v^2}{2}|\psi |^2 \end{aligned}$$

(34)

while the modified interaction Lagrangian is

$$\begin{aligned} {{\mathcal {L}}'}^I = {\mathcal {L}}^I + \frac{m_v^2}{2}|\psi |^2 = {\mathcal {L}}_3+ {\mathcal {L}}_4+ \frac{m_v^2}{2}|\psi |^2. \end{aligned}$$

The elements of the interacting algebra are now given in terms of two parameters $\lambda$ and $m_v$. More precisely, keeping $\mu$ fixed, as in (21), they are formal power series in $\sqrt{\lambda }$ with coefficients depending on $m_v^2$, which can be understood as a partial resummation of the original perturbative expansion. The advantage of this new expansion is in the fact that the coefficients remain finite in the adiabatic limit, when they are evaluated in the state representing the condensate at finite temperature. We furthermore observe that the principle of perturbative agreement discussed below implies that the final theory does not depend on this extra parameter $m_v$.

Let $H_{\psi ,\beta }$ be the symmetrized two-point function of the $\beta$-KMS given in (27), we observe that if $m_v$ is chosen to be sufficiently small, the interaction Lagrangian normal-ordered with respect to $H_{\psi ,\beta }$ is again convex. To see this in detail, let T be a time-ordering operator such that $TF = :F:_H$ where H is the distinguished Hadamard function constructed in (57). We have up to a choice of renormalization freedom [the lengthscale $\xi$ in (57) chosen in such a way that $:|\psi |^4:_H=:|\psi |^4:_{H_{\infty ,\psi }}$ where $H_{\infty ,\psi }$ is the symmetrized two-point function of the vacuum obtained taking the limit $\beta \rightarrow \infty$ in (27)]

$$\begin{aligned} T\left( \frac{1}{4}|\psi |^4\right)= & {} \frac{1}{4}:|\psi |^4:_{H_\beta }+ \frac{1}{2} (3 m_{\beta ,1}^2 + m_{\beta ,2}^2) :|\psi _1|^2:_{H_\beta } \nonumber \\&+ \frac{1}{2} (3 m_{\beta ,2}^2 + m_{\beta ,1}^2) :|\psi _2|^2:_{H_\beta } +C \end{aligned}$$

(35)

where C is a constant which can be discarded and the two thermal masses $m_{\beta ,i}$ have been computed above in (29) and (30)

$$\begin{aligned} m_{\beta ,1}^2 \doteq W(\psi _1^2), \quad m_{\beta ,2}^2 \doteq W(\psi _2^2). \end{aligned}$$

Hence $T(\frac{1}{4}|\psi |^4 -\frac{m_v^2}{2}|\psi |^2)$ remains convex, provided $m_v^2< \lambda (3 m_{\beta ,1}^2 + m_{\beta ,2}^2)$ and $m_v^2 < \lambda (3 m_{\beta ,2}^2 + m_{\beta ,1}^2)$. For this reason, it is expected that the stability properties of the theory are not altered adding the virtual masses $m_v$ in the free theory.

3.4 Condensate and Perturbative Agreement

We need to check that the Wick monomials in the interaction Lagrangian originally constructed over the linearized theory ${\mathcal {L}}_2$ are not corrected because of the new splitting. In other words, we prove that the principle of perturbative agreement holds also when a condensate is present. Let us recall the form of the equation of motion for ${{\tilde{\psi }}} = (\psi _1,\psi _2)$ given in (24)

$$\begin{aligned} D{{\tilde{\psi }}} = 0, \quad D = (\square -M^2) {\mathbb {I}} - \delta M^2 {{\varvec{\sigma }}}_3 - i 2\mu \partial _0 {{\varvec{\sigma }}}_2 \end{aligned}$$

and consider the preferred Hadamard function $H_{M^2,\delta M^2,\mu }$, with a lengthscale $\xi$, associated with this operator constructed in “Appendix C.”

We prove now that the time ordering operator $T_{M,\delta M,\mu }(F) = :F:_{H_{M^2,\delta M^2,\mu }}$ satisfies the principle of perturbative agreement. To this end, consider the $2\times 2$ matrix $\Psi =\{\psi _i\psi _j\}_{i,j\in \{1,2\}}$, following the discussion presented in “Appendix A” we want to prove that

$$\begin{aligned} \Delta \Psi = \gamma T_{0,0,\mu }\Psi - T_{M,\delta M,\mu } \Psi \end{aligned}$$

vanishes, where $\gamma$ is the map which intertwines $T_{0,0,\mu }(\psi _i(x)\psi _j(y))$ to $T_{M,\delta M,\mu }(\psi _i(x)\psi _j(y))$. Formally, indicating with the subscript c the quantities referred to the condensate $(M,\delta M,\mu )$ and with the subscript 0 those referred to the vacuum (0, 0, m) we have to compute

$$\begin{aligned} \Delta \Psi = \lim _{y\rightarrow x}\left( (H_F^0(x,y))_{\text {ren}} - H_F^c(x,y)\right) \end{aligned}$$

where $H_F^{c/0}$ are the time-ordered/Feynman propagator associated with the Hadamard functions $H^{c/0}$. By power counting, we notice that all the contributions larger than order two in $(D^c-D)$ are removed from $\Delta \Psi$ by renormalization. In order to check if there is a finite reminder after this renormalization, we analyze the form of the Hadamard singularity $H_{M^2,\delta M^2,\mu }$ given in “Appendix C”; we remove the contributions of order lower than the third in $x_i$ from $H_{M^2+x_1,\delta M^2+x_2,\mu +x_3}$ before computing the coinciding point limit. Let us recall the form of some Hadamard coefficient given in “Appendix C.” From Eqs. (59) and (60), we have

$$\begin{aligned} U =\cos (\mu x^0){\mathbb {I}} - i{{\varvec{\sigma }}}_2 \sin (\mu x^0) \end{aligned}$$

and

$$\begin{aligned} V_0 = -\frac{1}{2} U \left( (\mu ^2+M^2) {\mathbb {I}} + \delta M^2 \left( \frac{\sin (2\mu x^0)}{2\mu x^0}{{\varvec{\sigma }}}_3 +\frac{\cos (2\mu x^0)-1}{2\mu x^0} {{\varvec{\sigma }}}_1\right) \right) . \end{aligned}$$

Hence, in $H_{M^2+x_1,\delta M^2+x_2,\mu +x_3}$, $U/\sigma$ does not depend on $x_1$ and $x_2$ and the contributions in $x_3$ larger than the second order vanish in the coinciding point limits. Similarly, the contributions larger than second order vanish also in $V\log (\frac{\sigma }{\xi ^2})$. As for the mass perturbations, we thus have that the Wick monomials $\Psi ^n$ computed with respect to the Hadamard parametrix do not change under the action of the map which intertwines the time ordering constructed with two different sets of parameters $M,\delta M, \mu$.

3.5 Construction of the Condensate, Cluster Estimates

To construct the state at finite temperature over the condensate, we follow the construction given in [30] and summarized in Sect. 3.2. In particular, at fixed spatial cutoff h, the equilibrium state at inverse temperature $\beta$ can be constructed as in (11). Using the spatial translation invariance of the interacting hamiltonian and denoting by $\tau _{t,{\mathbf {x}}}$ the $*$-automorphisms realizing a spacetime translation of step $(t,{\mathbf {x}})$, we have, for any element A of ${\mathcal {A}}_I({\Sigma _\epsilon })$,

$$\begin{aligned} \omega ^{\beta ,V}_h(A) \doteq&\sum _{n}\int _{0\le u_1\le \cdots u_n\le \beta } d u_1\ldots d u_n\int _{{{\mathbb {R}}}^{3n}}\mathrm{d}^3{\mathbf {x}}_1\ldots \mathrm{d}^3{\mathbf {x}}_n h({\mathbf {x}}_1)\ldots h({\mathbf {x}}_n)\nonumber \\&\omega ^\beta _T\left( A;\tau _{iu_{1},{\mathbf {x}}_{1}}({\mathcal {H}}_I(0));\ldots ;\tau _{iu_n;{\mathbf {x}}_{n}}({\mathcal {H}}_I(0))\right) \end{aligned}$$

(36)

where $\omega ^\beta _T$ denotes the truncated n-point function of the state $\omega _{\beta ,\psi }$. Hence, in order to discuss the limit $h\rightarrow 1$ we need to control the decay for large spatial directions of the truncated n-point functions. We have actually the following theorem

Theorem 3.1

(Cluster expansions). Consider $A_i\in {\mathcal {A}}({\mathcal {O}})$ where ${\mathcal {O}}\subset B_R$ the open ball of radius R centered at the origin of the Minkowski spacetime and

$$\begin{aligned} F(u_1,{\mathbf {x}}_1;\ldots ;u_n,{\mathbf {x}}_n)\doteq \omega _T(A_0 ; \tau _{iu_1,{\mathbf {x}}_1}(A_1);\ldots ;\tau _{iu_n,{\mathbf {x}}_1}(A_n)). \end{aligned}$$

There exists a constant C such that

$$\begin{aligned} |F(u_1,{\mathbf {x}}_1;\ldots ;u_n,{\mathbf {x}}_n)| \le C e^{-\frac{m}{\sqrt{n}}r},\quad r=\sqrt{\sum _i |{\mathbf {x}}_i|^2} \end{aligned}$$

for $r>4cR$, uniformly in u for $0< u_1< \cdots< u_n <\beta$ with $\beta -u_n \ge \frac{\beta }{n+1}$.

Proof

Thanks to the decay property for large spatial separations of the locally smeared two-point functions given in Proposition D.2, the proof of this theorem can be done in a similar way as the proof of Theorem 3 of [30]. We recall here the main steps of that proof, and we adapt them to the case studied here.

The truncated n-point functions can be written as a sum over all possible connected graphs joining n points. We shall denote the set of connected graphs, without tadpoles, with $n+1$ vertices $V=\{0,\ldots , n\}$ as ${\mathcal {G}}_{n+1}^c$. Furthermore, for any $G\in {\mathcal {G}}_{n+1}$, E(G) denotes the set of edges of G. For any $l\in E(G)$, s(l) and r(l) denote the source and the range of l. A graph G is considered to be in ${\mathcal {G}}_{n+1}$ only if, for every l, $s(l)<r(l)$. Finally, $l_{ij}(G)$ is the number of lines connecting $i,j\in V$. With these definitions,

$$\begin{aligned} F(u_1,{\mathbf {x}}_1;\ldots ;u_n,{\mathbf {x}}_n)\doteq \sum _{G\in {\mathcal {G}}^c_{n+1}} \frac{1}{{\text {sym}(G)}} F_G(u_1,{\mathbf {x}}_1;\ldots ;u_n,{\mathbf {x}}_n) \end{aligned}$$

where $\text {sym}(G)=\prod _{i<j} l_{ij}(G)!$ is a numerical factor and

$$\begin{aligned}&F_G(u_1,{\mathbf {x}}_1;\ldots ;u_n,{\mathbf {x}}_n) \\&\quad \doteq \left( \prod _{0\le i<j\le n} \Gamma ^{ij} \right) \left. (A_0 \otimes \tau _{iu_1,{\mathbf {x}}_1}(A_1)\otimes \cdots \tau _{iu_n,{\mathbf {x}}_1}(A_n))\right| _{(\psi ^0,\ldots ,\psi ^n) = 0}. \end{aligned}$$

Furthermore,

$$\begin{aligned} \Gamma ^{ij} = \int \mathrm{d}^4x \mathrm{d}^4y \;{\mathcal {K}}(x-y) \frac{\delta }{\delta \psi ^i(x)} \otimes \frac{\delta }{\delta \psi ^j(y)} \end{aligned}$$

with the integral kernel ${\mathcal {K}}(x-y) = \omega _{\beta ,\psi }(\psi (x)\psi (y))$, given in terms of the thermal two-point function of the background theory (27). Furthermore, $\psi ^j=\psi ^j_1+i\psi ^j_2$ is the field configuration in the jth factor of the tensor product and the functional derivative $\frac{\delta }{\delta \psi ^j}$ acts on the jth factor of the tensor product. We have that

$$\begin{aligned} F_G(U,{\mathbf {X}}) \doteq \int dP \left( \prod _{l\in E(G)} e^{p_l^0(u_{s(l)}-u_{r(l)})} e^{i {\mathbf {p}}_l ({\mathbf {x}}_{s(l)}-{\mathbf {x}}_{s(l)})} \hat{{\mathcal {K}}}(p_l) \right) {{\hat{\Psi }}}(-P,P) \end{aligned}$$

where $U=(u_0,\ldots , u_n)$, ${\mathbf {X}}=({\mathbf {x}}_0, \ldots ,{\mathbf {x}}_n)$ with $u_0=0$ and ${\mathbf {x}}_0=0$, while $P=(p_1,\ldots ,p_{|E(G)|})$ and

$$\begin{aligned} \Psi (Z,Y) = \left. \left( \prod _{l\in E(G)} \frac{\delta }{\delta \psi ^{s(l)}(z_l)}\otimes \frac{\delta }{\delta \psi ^{r(l)}(y_l)} \right) (A_0 \otimes A_1\otimes \cdots A_n)\right| _{(\psi ^0,\ldots ,\psi ^n) = 0}. \end{aligned}$$

We observe that

$$\begin{aligned} \hat{{\mathcal {K}}}(p) = \left( \lambda _+(p) + \lambda _-(p)\right) (-\widehat{{\overline{D}}}) \end{aligned}$$

where $\lambda _+$ and $\lambda _-$ are the positive and negative frequency part

$$\begin{aligned} \lambda _+(p)&= \frac{1}{\omega _+^2-\omega _-^2} \left( \frac{\delta (p_0-\omega _+)}{2\omega _+}- \frac{\delta (p_0-\omega _-)}{2\omega _-}\right) \frac{1}{1-e^{-\beta p_0}}\\ \lambda _-(p)&= -\frac{1}{\omega _+^2-\omega _-^2} \left( \frac{\delta (p_0+\omega _+)}{2\omega _+}- \frac{\delta (p_0+\omega _-)}{2\omega _-}\right) \frac{1}{1-e^{-\beta p_0}}. \end{aligned}$$

Hence, separating the positive and negative contributions in $F_G$ we get

$$\begin{aligned} F_G(U,{\mathbf {X}})&= \sum _{P_2(E(G))}\int d{\mathbf {P}} \\&\quad \left( \prod _{l_+\in E_+(G)} e^{p_{l_+}^0(u_{s(l_+)}-u_{r(l_+)})} e^{i {\mathbf {p}}_{l_+} ({\mathbf {x}}_{s(l_+)}-{\mathbf {x}}_{r(l_+)})} \lambda _+(p_{l_+})(-\hat{{\overline{D}}}(p_{l_+}))\right) \\&\quad \cdot \left( \prod _{l_-\in E_-(G)} e^{p_{l_-}^0(u_{s(l_-)}-u_{r(l_-)})} e^{i {\mathbf {p}}_{l_-} ({\mathbf {x}}_{s(l_-)}-{\mathbf {x}}_{s(l_-)})} \lambda _-(p_{l_-})(-\hat{{\overline{D}}}(p_{l_-}))\right) \\&\quad {\hat{\Psi }}(-P,P) \end{aligned}$$

where the sum is taken over all possible partitions of E(G) in up to two sets $\{E_+(G),E_-(G)\}\in P_2(E(G))$. We proceed now splitting again these contributions over the two possible frequencies $\omega _{\pm }$. Hence, denoting by

$$\begin{aligned} \lambda _{++}(p)&\doteq \frac{1}{\omega _+^2-\omega _-^2} \left( \frac{\delta (p_0-\omega _+)}{2\omega _+}\right) \frac{1}{1-e^{-\beta p_0}}\\ \lambda _{+-}(p)&\doteq \frac{-1}{\omega _+^2-\omega _-^2} \left( \frac{\delta (p_0-\omega _-)}{2\omega _-}\right) \frac{1}{1-e^{-\beta p_0}}\\ \lambda _{-+}(p)&\doteq -\frac{1}{\omega _+^2-\omega _-^2} \left( \frac{\delta (p_0+\omega _+)}{2\omega _+}\right) \frac{1}{1-e^{-\beta p_0}}\\ \lambda _{--}(p)&\doteq \frac{1}{\omega _+^2-\omega _-^2} \left( \frac{\delta (p_0+\omega _-)}{2\omega _-}\right) \frac{1}{1-e^{-\beta p_0}} \end{aligned}$$

we have

$$\begin{aligned} F_G(U,{\mathbf {X}})&= \sum _{P_2(E(G))} \sum _{P_2(E_+(G))}\sum _{P_2(E_-(G))} \\&\quad \int d{P} \left( {\mathcal {Q}}_{++} \cdot {\mathcal {Q}}_{+-} \cdot {\mathcal {Q}}_{-+}\cdot {\mathcal {Q}}_{--}\right) {\hat{\Psi }}(-P,P) \end{aligned}$$

where

$$\begin{aligned}&{\mathcal {Q}}_{\sigma \sigma '} \doteq \left( \prod _{l\in E_{\sigma \sigma '}(G)} e^{p_{l}^0(u_{s(l)}-u_{r(l)})} e^{i {\mathbf {p}}_{l} ({\mathbf {x}}_{s(l)}-{\mathbf {x}}_{r(l)})} \lambda _{\sigma \sigma '}(p_{l})(-\hat{{\overline{D}}}(p_{l}))\right) \\&\quad \sigma ,\sigma '\in \{+,-\}. \end{aligned}$$

The function ${\hat{\Psi }}$ is an entirely analytic function which grows at most polynomially in every direction. We might thus integrate over all possible $p_0$ to get

$$\begin{aligned} F_G(U,{\mathbf {X}}) \doteq \sum _{P_2(E(G))} \sum _{P_2(E_+(G))}\sum _{P_2(E_-(G))}\int d{\mathbf {P}} \left( \tilde{{\mathcal {Q}}}_{++} \cdot \tilde{{\mathcal {Q}}}_{+-} \cdot \tilde{{\mathcal {Q}}}_{-+}\cdot \tilde{{\mathcal {Q}}}_{--}\right) \Phi ({\mathbf {P}}) \end{aligned}$$

where now

$$\begin{aligned} \Phi ({\mathbf {P}}) \doteq \left. {\hat{\Psi }}(-P,P)\right| _{p_0^{l_{\sigma \sigma '}} = \sigma \omega _{\sigma '}({\mathbf {p}}_l)} \end{aligned}$$

and

$$\begin{aligned}&\tilde{{\mathcal {Q}}}_{\sigma \sigma '} \\&\quad = (\sigma 1) (\sigma ' 1)\left( \prod _{l\in E_{\sigma \sigma '}(G)} \frac{e^{ \sigma \omega _{\sigma '}({\mathbf {p}}_{l})(u_{s(l)}-u_{r(l)})} e^{i {\mathbf {p}}_{l} ({\mathbf {x}}_{s(l)}-{\mathbf {x}}_{r(l)})} }{(\omega _+^2-\omega _-^2)2\omega _{\sigma '}} \frac{(-\hat{{\overline{D}}}(\sigma \omega _{\sigma '},{\mathbf {p}}_{l}))}{1-e^{-\sigma \beta \omega _{\sigma '}}} \right) \\&\quad = (\sigma ' 1) \left( \prod _{l\in E_{\sigma \sigma '}(G)} \frac{e^{-\left( \frac{(1-\sigma )}{2}\beta + \sigma (u_{r(l)}-u_{s(l)})\right) \omega _{\sigma '}} e^{i {\mathbf {p}}_{l} ({\mathbf {x}}_{s(l)}-{\mathbf {x}}_{r(l)})} }{(\omega _+^2-\omega _-^2)2\omega _{\sigma '}} \frac{(-\hat{{\overline{D}}}(\sigma \omega _{\sigma '},{\mathbf {p}}_{l}))}{1-e^{-\beta \omega _{\sigma '}}} \right) , \\&\qquad \sigma ,\sigma '\in \{+,-\}. \end{aligned}$$

Since, by hypothesis,

$$\begin{aligned} u_{i+1}>u_i,\quad \beta -u_{n} \ge \frac{\beta }{n+1} \end{aligned}$$

and $r(l)>s(l)$, we have that

$$\begin{aligned} e^{-\left( \beta - (u_{r(l)}-u_{s(l)})\right) \omega _{\sigma '}} \le e^{-\frac{n}{n+1}\beta \omega _{\sigma '}}. \end{aligned}$$

Hence,

$$\begin{aligned} {\tilde{\Phi }}({\mathbf {P}})\doteq \tilde{{\mathcal {Q}}}_{-+}\tilde{{\mathcal {Q}}}_{--}\Phi ({\mathbf {P}}) \end{aligned}$$

is rapidly decreasing, in every direction, because $F_G$ is a microcausal functional and $\Phi ({\mathbf {P}})$ is the restriction on a particular subdomain of ${\hat{\Psi }}(-P,P)$ which is an entire analytic function which grows at most polynomially. Hence, the negative frequencies are exponentially suppressed, and if directions containing only positive frequencies are considered, they are also rapidly decreasing by Proposition D.1. The integral over ${\mathbf {P}}$ can now be taken and we may apply Proposition D.2 to estimate the decay of the result of that integral. We obtain

$$\begin{aligned} |F_G(U,{\mathbf {X}})| \le c' \prod _{l\in E(G)} e^{- M_- \sqrt{|{\mathbf {x}}_{r(l)}-{\mathbf {x}}_{s(l)}|^2}}\le c' e^{- \frac{M_-}{\sqrt{n}} \sqrt{\sum _{i=1}^n |{\mathbf {x}}_{i}|^2}} \end{aligned}$$

where the constant $c'$ does not depend on $u_i$. In the last inequality, we used the fact that G is a connected graph, and thus, every $x_{i}$ can be reached from the origin $({\mathbf {x}}_0=0)$. Hence,

$$\begin{aligned} \sum _{l\in E(G)} \sqrt{{|{\mathbf {x}}_{r(l)}-{\mathbf {x}}_{s(l)}|^2}} \ge \text {max}_i{\sqrt{|{\mathbf {x}}_{i}|^2}} \ge \sqrt{\frac{1}{n}\sum _{i=1}^n |{\mathbf {x}}_{i}|^2} \end{aligned}$$

thus concluding the proof. $\square$

Theorem 3.2

Let $A \in {\mathcal {A}}_I({\mathcal {O}})$ where ${\mathcal {O}}\subset \Sigma _\epsilon$, the adiabatic limit

$$\begin{aligned} \begin{gathered} \omega ^{\beta ,V}(A)= \lim _{h\rightarrow 1}\sum _{n}\int _{0\le u_1\le \cdots u_n\le \beta } \mathrm{d} u_1\ldots \mathrm{d} u_n\int _{{{\mathbb {R}}}^{3n}}\mathrm{d}^3{\mathbf {x}}_1\ldots \mathrm{d}^3{\mathbf {x}}_n h({\mathbf {x}}_1)\ldots h({\mathbf {x}}_n) \\ \omega _T^\beta \left( A;\tau _{iu_{1},{\mathbf {x}}_{1}}(K);\ldots ; \tau _{iu_n;{\mathbf {x}}_{n}}(K)\right) , \end{gathered} \end{aligned}$$

where $K\doteq \lim _{h\rightarrow 1} {\mathcal {H}}_I(0)$, exists in the sense of perturbation theory and defines an equilibrium state for the interacting theory.

Proof

Since ${\mathcal {O}}$ is of compact support, it exists and $R>0$ such that the open ball $B_R$ centered in the origin of Minkowski spacetime contains ${\mathcal {O}}$, namely ${\mathcal {O}}\subset B_R$. Furthermore, thanks to the temporal cutoff $\chi$, and in view of the causal properties of the Bogoliubov map, $K = \lim _{h\rightarrow 1} {\mathcal {H}}_I$ is supported in $B_R$ for a sufficiently large R.

Consider the nth order contribution in the sum defining $\omega ^{\beta ,V}_h$ given in (36)

$$\begin{aligned} \begin{aligned} \Omega _{n,h}(A)&\doteq \int _{0\le u_1\le \cdots u_n\le \beta } d u_1\ldots d u_n \\&\quad \int _{{{\mathbb {R}}}^{3n}}\mathrm{d}^3{\mathbf {x}}_1\ldots \mathrm{d}^3{\mathbf {x}}_n h({\mathbf {x}}_1)\ldots h({\mathbf {x}}_n) F(u_1,{\mathbf {x}}_1;\ldots ;u_n,{\mathbf {x}}_n) \end{aligned} \end{aligned}$$

(37)

where

$$\begin{aligned} F(u_1,{\mathbf {x}}_1;\ldots ;u_n,{\mathbf {x}}_n)\doteq \omega ^\beta _T\left( A;\tau _{iu_{1},{\mathbf {x}}_{1}}(K);\ldots ;\tau _{iu_n;{\mathbf {x}}_{n}}(K))\right) , \quad A\in {\mathcal {A}}_I({\mathcal {O}}). \end{aligned}$$

To apply the results of Theorem 3.1, we observe that if R is sufficiently large ${\mathcal {H}}_I(0) \in {\mathcal {A}}_I(B_R)$, furthermore, the form of the integration domain of the u variables as given in (37) is such that

$$\begin{aligned} 0 \le u_1 \le \cdots \le u_n \le \beta . \end{aligned}$$

(38)

Using the KMS condition, we might restrict attention to the case where $\beta - {u_{n}} \ge \frac{\beta }{n+1}$. In fact, if this is not the case, there must exist an m for which $u_{m}-u_{m-1} \ge \frac{\beta }{n+1}$. Actually, for $A_i\in {\mathcal {A}}_I({\mathcal {O}})$ by the KMS condition we have that

$$\begin{aligned}&\omega _T^\beta (\tau _{iu_0}(A_0) ; \tau _{iu_1,{{\mathbf {x}}}_1}(A_1); \ldots ; \tau _{iu_n,{{\mathbf {x}}}_n}(A_n))\\&\quad =\omega _T^{\beta }(\tau _{iu_{m},{\mathbf {x}}_{m} } (A_{m}) ; \ldots ; \tau _{iu_n,{{\mathbf {x}}}_n}(A_n) \\&\qquad \otimes \tau _{i\beta +iu_0}(A_0) ; \ldots ; \tau _{i\beta +iu_{m-1},{{\mathbf {x}}}_{m-1}}(A_{m-1}) ); \end{aligned}$$

hence, we might now consider

$$\begin{aligned}&F'(v_1,{\mathbf {y}}_1;\ldots ;v_n,{\mathbf {y}}_n) \\&\quad \doteq \omega ^\beta _T(K ; \tau _{iv_1,{\mathbf {y}}_1}(K);\ldots ; \tau _{iv_{n-m},{\mathbf {y}}_{n-m}}(K); \tau _{iv_{n-m+1},{\mathbf {y}}_{n-m+1}}(A_0); \\&\qquad \tau _{iv_{n-m+2},{\mathbf {y}}_{n-m+2}}(K);\ldots \tau _{iv_n,{\mathbf {y}}_v}(K)) \end{aligned}$$

in place of F. In fact, the previous equality obtained with the KMS condition together with translation invariance of the state implies that

$$\begin{aligned} F(u_1,{\mathbf {x}}_1;\ldots ;u_n,{\mathbf {x}}_n) = F'(v_1,{\mathbf {y}}_1;\ldots ;v_n,{\mathbf {y}}_n) \end{aligned}$$

if

$$\begin{aligned} ({\mathbf {y}}_1,\ldots ,{\mathbf {y}}_n) = ({\mathbf {x}}_{m+1}-{\mathbf {x}}_m,\ldots , {\mathbf {x}}_{n}-{\mathbf {x}}_m, -{\mathbf {x}}_m,{\mathbf {x}}_1-{\mathbf {x}}_m,\ldots ,{\mathbf {x}}_{m-1}-{\mathbf {x}}_m) \end{aligned}$$

and

$$\begin{aligned} (v_1,\ldots ,v_n) = (u_{m+1}-u_m,\ldots , u_{n}-u_m,\beta -u_m,\beta +u_1-u_m,\ldots ,). \end{aligned}$$

The arguments of the function $F'$ have the desired property, actually $\beta - v_n= u_m - u_{m-1} \ge \beta /(n+1)$. We might thus use $F'$ in place of F, because the integration over the u variables is over a compact set and because the points where $u_i=u_j$ for some $i\ne j$ form a zero measure set. Hence, Theorem 3.1 implies that the integral over ${\mathbf {x}}_i$ can be taken for all i to conclude the proof. $\square$

4 Spontaneous Symmetry Breaking and the Goldstone Theorem

The model we are considering possesses an internal U(1) symmetry. Actually, the Lagrangian is invariant under transformations

$$\begin{aligned} \varphi = U(\theta )\varphi \doteq e^{i\theta } \varphi \end{aligned}$$

where $\theta \in [0,2\pi ]$. However, in the state which describes the condensate, this symmetry is spontaneously broken. From the Goldstone theorem, we expect that a massless (gapless) mode is present in the model. This observation is in contrast with the analysis discussed in the previous section. Actually there, all the fields in the linearized theory were assumed to be massive. Notice that once the background $\phi$ is fixed in the decomposition $\varphi =e^{-i\mu x^0}(\phi +\psi )$, the Lagrangian for the linearized theory is not invariant under U(1) transformations

$$\begin{aligned} \psi \rightarrow e^{i\theta } \psi + (e^{i\theta }-1)\phi ; \end{aligned}$$

hence, the fact that both linearized fields $\psi _i$ are massive is not in contrast with the Goldstone theorem. Furthermore, if the Goldstone theorem holds for the full theory, this would imply that at least one gapless mode should exist if the full perturbation series is considered.

In the case of thermal theories, the proof of Goldstone theorem is not completely straightforward as for theories at zero temperature because the original proof makes use of Lorentz invariance [34] (see also the work of Jona-Lasinio using effective action methods in [45]). The equilibrium states are, however, not Lorentz invariant because of the presence of a preferred time direction in the KMS condition. Furthermore, even if a gapless mode exists, the particle content of the gapless mode is not immediately evident, as discussed by Bros and Buchholz in [11].

The presence of Goldstone modes at finite temperature has been discussed in [48] using effective action methods. Based on the analysis of Swieca [69], a proof of the Goldstone theorem without using Lorentz invariance has been given by Morchio and Strocchi in [54], see also the book [65] for the application of similar ideas for the analysis of the case of finite temperature. Furthermore, the analysis of the slow decay of large spatially separated correlation functions in the presence of spontaneous symmetry breaking is discussed in [44]. However, when a nontrivial background is present as for the case of Bose–Einstein condensation, we don’t expect that the presence of a gapless mode is directly related to the clustering properties of the correlation functions for large spatial separation. As an example, consider the two-point function of the state $\omega ^\pm _{\beta ,c}$ discussed in Sect. 3.1 in the limit of vanishing temperature, namely $\beta \rightarrow \infty$. The obtained state is the composition of the massive vacuum $\omega _0$ with the map $\gamma ^\pm _c$ given in (19). Even if one of the modes in the two-point function of $\omega ^\pm _{\beta ,c}$ is gapless, the clustering properties of $\omega ^\pm _{\beta ,c}$ are equivalent to the one of the vacuum because $\gamma ^\pm _c$ does not change the localization of the observables.

The mentioned proofs cannot be directly applied for perturbatively constructed theories; for this reason in the next section we shall give a proof of the validity of Goldstone theorem which can hold in our setting. We shall actually follow Swieca’s proof without making use of Lorentz invariance.

4.1 Proof of the Goldstone Theorem

Here we would like to give a proof of the validity of Goldstone theorem at finite temperature in the presence of a condensate. For this purpose, we observe that the U(1) invariance of the Lagrangian density ${\mathcal {L}}$ for $\varphi =\varphi _1+i\varphi _2$ is such that

$$\begin{aligned} (\varphi _1,\varphi _2)\rightarrow R(\theta )(\varphi _1,\varphi _2) \end{aligned}$$

(39)

where R is a rotation of an angle $\theta$. Its infinitesimal version is

$$\begin{aligned} \varphi _m \rightarrow \varphi _m+\epsilon t_{mn} \varphi _n \end{aligned}$$

where t is the antisymmetric metric

$$\begin{aligned} t= \begin{pmatrix} 0&{}\quad -1\\ 1&{}\quad 0 \end{pmatrix}. \end{aligned}$$

Associated with the symmetry which is spontaneously broken in the state $\omega$, there is a current J which is conserved. This current is defined as

$$\begin{aligned} J^\mu \doteq \frac{\delta {\mathcal {L}}}{\delta \partial _\mu \varphi _m }t_{mn}\varphi _n = i ( {\overline{\varphi }}\partial ^\mu {\varphi } -\partial ^\mu {\overline{\varphi }}\varphi ). \end{aligned}$$

(40)

By Noether theorem, the action possesses the desired U(1)-symmetry if and only if the current J is conserved, namely if $\nabla _\mu J^\mu =0$. Following [46], we can now introduce a regularized charge operator associated with the current density $J^0$ introduced above. Let $f\in C^{\infty }_0({\mathbb {R}})$ be a time cutoff with $\text {supp}{f}\in (-\epsilon ,\epsilon )$, $f\ge 0$ and $\Vert f\Vert _1=1$. Furthermore, $g\in C^\infty _{0}({\mathbb {R}}^3)$ is a space cutoff, $g({\mathbf {x}})=1$ for ${\mathbf {x}}<1$. The regularized charge operator associated with J can be seen as the large R limit of

$$\begin{aligned} Q_R \doteq \int \mathrm{d}^4x f(x_0) g\left( \frac{{\mathbf {x}}}{R}\right) J^0(x). \end{aligned}$$

(41)

The charge operator can be used to implement the infinitesimal U(1) transformation of the field⁴

$$\begin{aligned} \lim _{R\rightarrow \infty }\left[ Q_R,\varphi _m(t,{\mathbf {y}})\right] = t_{mn} \varphi _n(t,{\mathbf {y}}). \end{aligned}$$

(42)

Hence, in a state where the symmetry is spontaneously broken, namely when $\omega (\varphi _n(0))=\phi _n\ne 0$ for some n,

$$\begin{aligned} \lim _{R\rightarrow \infty } \omega ([ Q_R,\varphi _n(0)]) =t_{nm} \omega ({\varphi }_{m}(0)) =t_{nm}{{\phi }}_{m}. \end{aligned}$$

(43)

Notice that, in view of the support properties of f, g, and of the conservation of the current J, for R sufficiently large we have that $\omega ([ Q_R,\varphi _n(0)])$ is constant. Hence, the limit $R\rightarrow \infty$ can be safely taken and the final result does not depend on the particular form of g.

We are now ready to state the Goldstone theorem in the following form

Theorem 4.1

(Goldstone). Consider a complex scalar quantum field $\varphi$ whose Lagrangian density ${\mathcal {L}}$ possesses an U(1)-symmetry generated by the current J given in (40). Consider the distribution

$$\begin{aligned} G_n(x)\doteq \omega ([J^0(x),\varphi _n(0)]). \end{aligned}$$

If the symmetry generated by J is spontaneously broken in the state $\omega$, namely $\omega (\varphi _n(0))=\phi _n\ne 0$ for some n, then in the spectrum of $G_n$ there is a zero frequency (gapless) contribution at vanishing momentum; namely ${\hat{G}}_n$, the Fourier transform of $G_n$, is such that

$$\begin{aligned} \lim _{{\mathbf {p}}\rightarrow 0} {\hat{G}}_n(p_0,{\mathbf {p}}) = \delta (p_0) t_{nm}\phi _m \end{aligned}$$

(44)

in the sense of distribution.

Proof

The support properties of the distribution $F^i_n(x) \doteq \omega ([J^i(x),\varphi _n(0)])$ and invariance under spatial rotations imply that there exists a distribution ${\mathcal {F}}_n$ such that

$$\begin{aligned}&\int \mathrm{d}^4x f(x_0) g\left( \frac{{\mathbf {x}}}{R}\right) \omega ([J^i(x),\varphi _n(0)]) \\&\quad =\int {d}^4p {\hat{f}}(p_0) {\hat{g}}({\mathbf {p}}R)R^3 {\mathbf {p}}^i {\mathcal {F}}_n(p^0,|{\mathbf {p}}|)\\&\quad =\int {d}^4p {\hat{f}}(p_0) {\hat{g}}({\mathbf {p}}) \frac{{\mathbf {p}}^i}{R} {\mathcal {F}}_n\left( p^0,\frac{|{\mathbf {p}}|}{R}\right) , \quad i\in \{1,2,3\} \end{aligned}$$

where f and g are chosen as in the definition of regularized charge $Q_R$ given in (41). Causality implies that for R sufficiently large, the left-hand side does not depend on R. Hence, $\int dp_0 {\hat{f}}(p_0) {{\mathbf {p}}^i} {\mathcal {F}}_n(p^0,{|{\mathbf {p}}|})$ must be bounded near ${\mathbf {p}}=0$ and as a consequence of this fact it holds that

$$\begin{aligned}&\lim _{R\rightarrow \infty }\int \mathrm{d}^4x f(x_0) g\left( \frac{{\mathbf {x}}}{R}\right) \sum _{i=1}^3\nabla _i\omega ([J^i(x),\varphi _n(0)]) \\&\quad = \lim _{R\rightarrow \infty } \int {d}^4p {\hat{f}}(p_0) {\hat{g}}({\mathbf {p}}) \frac{|{\mathbf {p}}|^2}{R^2} {\mathcal {F}}_n\left( p^0,\frac{|{\mathbf {p}}|}{R}\right) =0. \end{aligned}$$

Hence, current conservation furnishes a condition for $G_n$, namely

$$\begin{aligned}&\lim _{R\rightarrow \infty }\int \mathrm{d}^4x f(x_0) g\left( \frac{{\mathbf {x}}}{R}\right) \sum _{\nu =0}^3\nabla _\nu \omega ([J^\nu (x),\varphi _n(0)]) \\&\quad = \lim _{R\rightarrow \infty }\int \mathrm{d}^4x f(x_0) g\left( \frac{{\mathbf {x}}}{R}\right) \nabla _0\omega ([J^0(x),\varphi _n(0)]) \\&\quad = \lim _{R\rightarrow \infty }\int \mathrm{d}^4p {\hat{f}}(p_0) {\hat{g}}({\mathbf {p}}) p_0 {\hat{G}}_n\left( p_0,\frac{{\mathbf {p}}}{R}\right) \mathrm{d}^4p =0. \end{aligned}$$

Another constraint on the form of $G_n$ is given by (43), namely

$$\begin{aligned} \lim _{R\rightarrow \infty }\int \mathrm{d}^4p {\hat{f}}(p_0) {\hat{g}}({\mathbf {p}}) {\hat{G}}_n\left( p_0,\frac{{\mathbf {p}}}{R}\right) \mathrm{d}^4p = t_{nm} \phi _m {\hat{f}}(0). \end{aligned}$$

Both conditions imply that (44) holds in the sense of distributions. $\square$

As discussed in [11], no direct particle interpretation can be inferred from (44). Actually, the singularity in ${\hat{G}}_n$ can be proved to exist only at ${\mathbf {p}}=0$ and not on the whole null cone. For a particle interpretation of this fact, we refer to the paper [11].

4.2 Analysis of the Validity of Goldstone Theorem in Perturbation Theory

We have seen that Goldstone theorem in the form stated in Theorem 4.1 holds for a quantum scalar field theory if the corresponding Lagrangian density ${\mathcal {L}}$ is invariant under the U(1) transformations given in (39), and if this symmetry is spontaneously broken in a state $\omega$. Hence, in a generic quantum field theory, we can apply Goldstone theorem if the current J given in (40) is conserved and if $\omega (\varphi _n(0))=\phi _n\ne 0$ for some n.

We now check that at linear order some of the desired hypotheses are not satisfied. In particular, we see that in the linearized theory the internal U(1) symmetry is explicitly broken. Actually, the current J defined as in (40), whose time component has the explicit form (31), is not conserved. To check this fact, notice that

$$\begin{aligned} \begin{aligned} \partial ^\mu J_\mu&= 2(\psi _2\square \psi _1-\psi _1\square \psi _2) -2 \phi \square \psi _2 -4 \mu (\phi {\dot{\psi }}_1+\psi _1{\dot{\psi }}_1 +\psi _2{\dot{\psi }}_2)\\&= 2\psi _1\psi _2(M_1^2-M_2^2) - 2\phi M_2^2 \psi _2 \end{aligned} \end{aligned}$$

(45)

where we have used the equation of motion given at linear order (22). If now $\mu ^2>m^2$, we have that $\phi \ne 0$ because $\lambda \phi ^2 = (\mu ^2-m^2)$. Furthermore, $M_1^2=(m^2-\mu ^2) + 3\lambda \phi ^2$, and $M_2^2=(m^2-\mu ^2) + \lambda \phi ^2$, and thus, even if $M_2=0$, we have that $M_1^2-M_2^2\ne 0$, and hence, $\partial ^\mu J_\mu \ne 0$.

We pass now to analyze the interacting case. We observe that if the full interacting equation of motion is used in evaluating $\partial ^\mu J_\mu$, namely taking into account ${\mathcal {L}}_3$ and ${\mathcal {L}}_4$, we have that equation (45) needs to be changed to

$$\begin{aligned} \partial ^\mu J_\mu = 2\psi _1\psi _2(M_1^2-M_2^2-2\lambda \phi ^2) - 2\phi M_2^2 \psi _2. \end{aligned}$$

Now both $M_2=0$ and $M_1^2-M_2^2-2\lambda \phi ^2=0$, and hence, the symmetry is not explicitly broken in the full classical theory. Notice that this analysis does not depend on the splitting between free and interacting Lagrangian; hence, the eventual thermal mass contributions do not alter this analysis.

We now check that conservation of the current J holds in the case of interacting quantum field theory treated with perturbation methods. The very same analysis holds up to the quotient with respect to the free equation of motion. Actually, the Schwinger–Dyson equation implies that

$$\begin{aligned} R_V\left( \frac{\delta S}{\delta \varphi } \right) = \frac{\delta S_0}{\delta \varphi } \ne 0 \end{aligned}$$

where S is the action constructed with the full Lagrangian density ${\mathcal {L}}$ and $S_0$ is the action of the linearized theory, constructed with the Lagrangian density ${\mathcal {L}}_2$. If expectation values on a quantum state are considered, we have

$$\begin{aligned} \omega \left( R_V\left( \frac{\delta S}{\delta \varphi } \right) \right) = 0. \end{aligned}$$

We now analyze in the same spirit the conservation of the current density J. To this end, we observe that

$$\begin{aligned} R_V(\partial ^\mu J_\mu )= & {} i R_V({\overline{\varphi }}\square \varphi - \square {\overline{\varphi }}\varphi ) = R_V\left( {{\overline{\varphi }}}\frac{\delta S}{\delta {{\overline{\varphi }}}}- \varphi \frac{\delta S}{\delta \varphi } \right) \nonumber \\= & {} R_V\left( {{\overline{\varphi }}} \cdot _T \frac{\delta S}{\delta {{\overline{\varphi }}}}- \varphi \cdot _T \frac{\delta S}{\delta \varphi } \right) \end{aligned}$$

(46)

where the last equality holds because of the properties of the time-ordered product and ${\overline{S}}=S$. We stress that the divergences present in ${{\overline{\varphi }}}(x)\cdot _T \delta S/\delta {{\overline{\varphi }}}(x)$ and proportional to $\lim _{y\rightarrow x}\Delta _F(x,y) \frac{\delta ^2 S}{\delta \varphi (y)\delta {{\overline{\varphi }}}(x)}$ are canceled by the subtraction of $\varphi \cdot _T \delta S/\delta \varphi$. To proceed, we need the Master Ward Identity which can be shown to hold for this theory without anomalies. Master Ward Identities have been discussed in [26, 27, 39]. For our purposes, we start with equation (70) of [31]. Rewriting it in our context, we obtain the desired equation⁵

$$\begin{aligned} R_V\left( {{\overline{\varphi }}}\cdot _T \frac{\delta S}{\delta {{\overline{\varphi }}}}- \varphi \cdot _T\frac{\delta S}{\delta \varphi } \right) = R_V({{\overline{\varphi }}}) \star \frac{\delta S_0}{\delta {{\overline{\varphi }}}}- R_V(\varphi )\star \frac{\delta S_0}{\delta \varphi }. \end{aligned}$$

(47)

Notice that the divergences on the right-hand side of the previous equality due to the pointwise multiplication with the quantum product vanishes because $\delta S_0/\delta \varphi$ is in the ideal of the linear equation of motion. Hence, in any quantum state

$$\begin{aligned} \omega (R_V(\partial ^\mu J_\mu )) = \omega \left( R_V({{\overline{\varphi }}}) \star \frac{\delta S_0}{\delta {{\overline{\varphi }}}}- R_V(\varphi )\star \frac{\delta S_0}{\delta \varphi } \right) = 0 \end{aligned}$$

We conclude that the current J constructed with the interaction quantum scalar field is conserved in the sense of perturbation theory up to the ideal describing the equation of motion. Finally, we observe that even if this last observation has been given in terms of the field $\varphi$, since relations (46) and (47) are algebraic relations, the very same analysis holds for the fluctuations $\psi _1,\psi _2$.

We also observe that on the state we are considering $\omega (\varphi _1) \ne 0$ in the sense of perturbation theory. Actually, $\phi = \frac{1}{\sqrt{\lambda }} \sqrt{\mu ^2-m^2}$, and thus, ${\mathcal {L}}_3$ is of order $\sqrt{\lambda }$ while ${\mathcal {L}}_4$ is of order $\lambda$ ($m_v^2$ needs to be chosen smaller than $\lambda c$). This implies that $\omega (R_V \psi _1)$ is at least of order O(1) in $\lambda$, and hence, it cannot totally cancel/compensate $\phi$.

We thus have that the hypotheses of Theorem 4.1 hold in the sense of perturbation theory. Furthermore, we also notice that all the identities in the proof hold in the sense of perturbation theory. Hence, we conclude that the thesis of that theorem holds in the sense of perturbation theory.

Appendices

The Principle of Perturbative Agreement

Thanks to the principle of perturbative agreement [24, 42] the theories obtained with the two different splittings given in (12) and (13) of the Lagrangian density

$$\begin{aligned} {\mathcal {L}} = -\frac{1}{2}\partial \phi \partial \phi - \frac{\lambda }{4} \phi ^4 \end{aligned}$$

are equivalent. In particular, the principle of perturbative agreement holds because it exists a map $\gamma$ which intertwines the time-ordered products constructed in theories with different masses m and $m'$ which will be later set to 0. We shall here denote by $T_m$ the time ordering operator of a massive theory $\square -m^2$. Hence,

$$\begin{aligned} \gamma : T_m {\mathcal {F}}_{\text {mloc}} \rightarrow T_{m'} {\mathcal {F}}_{\text {mloc}} \end{aligned}$$

where here ${\mathcal {F}}_{\text {mloc}}$ is the commutative algebra of multilocal functionals (with respect to the pointwise (classical) product). The nontrivial nature of $\gamma$ can be seen from the fact that the time-ordered products $T_m$ and $T_{m'}$ are different. At the same time, the perturbative agreement requires that the common part of ${\mathcal {L}}_I$ and ${\mathcal {L}}'_I$, or more generally all Wick powers, is left invariant up to ordinary renormalization freedom by the map $\gamma$. We shall here see that this is the case provided the time-ordering operator $T_m$ is constructed with the preferred symmetrized Hadamard singularity $H_m$ of $\square -m^2$ as in (2) with $W=0$ for some fixed lengthscale $\xi$, namely $T_m(F)=:F:_{H_m}$.

For completeness, we discuss this invariance under $\gamma$ in the simple case of the Wick square

$$\begin{aligned} \Phi ^2(f) = \int \phi ^2(x) f(x) \mathrm{d}^4x. \end{aligned}$$

We want to compare $T_{m'} \Phi ^2(f)$ with $T_m \Phi ^2(f)$. To this end, we observe that

$$\begin{aligned} T_{m'}\Phi ^2(x) = \lim _{y\rightarrow x}T_{m'}(\varphi (x)\varphi (y)) - \Delta _{F,{m'}}(x,y), \end{aligned}$$

and we study

$$\begin{aligned} \Delta \Phi ^2=\gamma T_{m'}\Phi ^2 - T_{m}\Phi ^2. \end{aligned}$$

Suppose that, on regular functionals, $T_{m'}$ and $T_{m}$ are constructed starting with the Feynman propagators $\Delta _{F,{m'}}(x,y)$ and $\Delta _{F,{m}}(x,y)$; hence, we have

$$\begin{aligned} \Delta \Phi ^2(x)&=\lim _{y\rightarrow x}\left( \gamma T_{m'}(\varphi (x)\varphi (y)) -T_{m}(\varphi (x)\varphi (y)) - \Delta _{F,{m'}}(x,y) + \Delta _{F,{m}}(x,y)\right) . \end{aligned}$$

Furthermore, the principle of perturbative agreement (see, e.g., [24] Theorem 3.2 where $\gamma$ is denoted by $\beta$) implies that $\gamma \varphi = \varphi$ and that

$$\begin{aligned} \gamma T_{m'}(\varphi (x)\varphi (y)) = T_{m}(\varphi (x)\varphi (y)), \end{aligned}$$

and thus,

$$\begin{aligned} \Delta \Phi ^2(x)&=\lim _{y\rightarrow x}\left( - \Delta _{F,{m'}}(x,y) + \Delta _{F,{m}}(x,y)\right) . \end{aligned}$$

It remains to compare $\Delta _{F,{m'}}$ with $\Delta _{F,m}$. Considering the form of the Hadamard singularities, we observe that

$$\begin{aligned} \lim _{y\rightarrow x}\left( \Delta _{F,{m'}}(x,y) - \Delta _{F,{m}}(x,y) \right) = \lim _{y\rightarrow x} ( V_{m'} - V_{m})\log \left( \frac{\sigma _\epsilon }{\mu ^2}\right) + W(x,y) \end{aligned}$$

where $V_m$ is the standard Hadamard V coefficient, $\sigma _\epsilon$ is the regularized squared geodesic distance, $\mu$ is a length scale and W is a smooth reminder. In the coinciding point limit, $( V_{m'} - V_{m})$ is proportional to $\delta m^2 = {m'}^2-{m}^2$; hence, the limits present in the previous equation are divergent. However, the kind of divergences present in that coinciding point limit is not different than the standard divergences present in $R_Q(\varphi ^2)$ with

$$\begin{aligned} Q \doteq \int g \delta m^2 \varphi ^2 \mathrm{d}^4x, \end{aligned}$$

if the time-ordered product of $T_m(Q\varphi ^2)$ is constructed without making use of a correct renormalization prescription. We may thus avoid (renormalize) them in a similar way as the divergences present in the naive construction of $T_m$ are resolved. To this end, we notice that if we consider the Feynman operator as operator on functions, we have that

$$\begin{aligned} \Delta _{F,m'} = \lim _{g\rightarrow 1}\sum _{n\ge 0}\Delta _{F,m} \left( -g \delta m^2 \Delta _{F,m}\right) ^n. \end{aligned}$$

We observe, by power counting, that the divergent contributions in the sum are the contributions $n=0$ and $n=1$. The contributions $n=0$ are removed by the subtraction of $\Delta _{F,m}$ while the divergences present in the contribution $n=1$ are similar to the divergences present in $T_m(QQ)$, and thus, they can be treated with renormalization theory. Actually, at order 1 in $\delta m^2$ the remaining contribution $\Delta \Phi ^2_{(1)}$ in the difference $\gamma T_{m'}\Phi ^2-T_{m}\Phi ^2$ is thus

$$\begin{aligned} \Delta \Phi ^2_{(1)} = -\lim _{g\rightarrow 1}\int (\Delta _{F,m}^2(x-y) )_{\text {ren}} \delta m^2 g(y)\mathrm{d}^4y. \end{aligned}$$

We regularize the product of Feynman propagators in the following way

$$\begin{aligned}&(\Delta _{F,m}^2)_{\text {ren}}(x) = (\square + a )\int _{(2m)^2}^\infty dM^2 \frac{\rho _2}{M^2+a} i\Delta _{F,M}(x), \nonumber \\&\rho _2=\frac{1}{16\pi ^2} \sqrt{1-\frac{4m^2}{M^2}} \end{aligned}$$

(48)

where a is a parameter which takes into account the known renormalization freedom. If we fix it to 0, we have that

$$\begin{aligned} -\Delta \Phi ^2_{(1)}&= \lim _{g\rightarrow 1}\frac{\delta m^2}{(2\pi )^{-4}}\int { d}^4p \;{\hat{g}}(p) ({\hat{\Delta }}_{F,m}^2)_{\text {ren}}(p) \\&=\lim _{g\rightarrow 1} \frac{\delta m^2}{(2\pi )^{-4}} \int {d}^4p \;{\hat{g}}(p) (-p^2)\int _{(2m)^2}^\infty dM^2 \frac{\rho _2}{M^2} \frac{1}{-p^2+M^2+i\epsilon } =0 \end{aligned}$$

where in the last step we used the fact that ${\hat{g}}$ tends to $\delta (p)$ in the limit $g\rightarrow 1$. Higher-order contributions can be computed directly, and they give a finite result.

Instead of performing these explicit computations, for our purposes, we need to understand the origin of these contributions with the following observation. Actually, the essential steps in the previous computation are the following: We have computed the expansion in powers of $\delta m^2$ of $\Delta _{F,m'}(x-y)$, we have removed the contributions of order 0 and 1, and we have eventually taken the coinciding point limits. In formulae

$$\begin{aligned} -\Delta \Phi ^2 = \lim _{y\rightarrow x} \left( \Delta _{F,m'}(x,y)-\Delta _{F,m}(x,y)- \left. \frac{\partial }{\partial \delta m^2}\Delta _{F,m'}(x,y)\right| _{\delta m^2=0}\delta m^2 \right) .\nonumber \\ \end{aligned}$$

(49)

We prove now that if one starts with the time-ordered propagators constructed with the preferred Hadamard functions H instead of the Minkowski vacuum $\Delta \Phi ^2$ vanishes, we have actually that

$$\begin{aligned} H_{F,m'} = \frac{U}{\sigma _\epsilon } + V\log \left( \frac{\sigma _\epsilon }{\mu ^2}\right) \end{aligned}$$

and expanding V in powers of $\sigma$ we have

$$\begin{aligned} V= c \frac{{m'}^2}{m'\sqrt{\sigma }}I_1(m'\sqrt{\sigma }) = \frac{c}{2}\left( m'^2 + \frac{1}{8}{m'}^4 \sigma + \cdots \right) . \end{aligned}$$

Hence, we conclude that

$$\begin{aligned} -\Delta \Phi ^2 = \lim _{y\rightarrow x} \left( H_{F,m'}(x,y)-H_{F,m}(x,y)- \left. \frac{\partial }{\partial \delta m^2}H_{F,m'}(x,y)\right| _{\delta m^2=0} \delta m^2\right) = 0. \end{aligned}$$

The same argument essentially holds also on any curved background, so in the case of mass renormalization we see that there is a choice of renormalization freedom [$a=0$ in (48)] such that

$$\begin{aligned} \gamma T_{m'}\Phi ^2 =T_m\Phi ^2. \end{aligned}$$

where now the Wick powers $T_m\Phi ^2 = :\Phi ^2:_{H_m}$ are constructed (regularized) with respect to the distinguished Hadamard singularity $H_m$ and not with respect to the symmetric part the two-point function of a state. Finally, we observe that the same results hold also when Wick monomial of higher order is considered, namely $\gamma T_{m'}\Phi ^n=T_{m}\Phi ^n$.

Nonrelativistic Limit of the Free Complex Scalar Field

In Sect. 3.1, we have analyzed the possible KMS states with $\beta >0$ and chemical potential $\mu$. In particular, for $\beta >0$ and chemical potential $|\mu |<m$ we have $\omega _{\beta ,\mu }$, for $\beta >0$ and chemical potential $\mu =\pm$ and a condensate c we have $\omega _{\beta ,c}^{\pm }$. We discuss in this appendix the nonrelativistic limit of these states. We shall, furthermore, see that the charge density converges to the particle density under that limit.

The nonrelativistic limit is obtained for temperatures $T=\beta ^{-1}$ which are small compared to m and for velocities $v\ll 1$ such that $mv^2/T=O(1)$. The chemical potential is near to m, $(m-\mu )/T=O(1)$. We set ${\bar{\mu }}=\lambda ^{-2}(m-\mu )$ and

$$\begin{aligned} \psi _{\lambda }(t,{\mathbf {x}})=\sqrt{2m\lambda ^3}\ \varphi (\lambda ^2 t,\lambda {\mathbf {x}})e^{-it\lambda ^2m} \end{aligned}$$

(50)

and find for the first class of KMS states

$$\begin{aligned}&\omega _{\lambda ^2\beta ,m-\lambda ^{-2}\bar{\mu }}(\psi _{\lambda }^*(t,x)\psi _{\lambda }(t',x')) \nonumber \\&\quad =\frac{1}{(2\pi )^3}\int {d}^4 p\frac{\delta (p^2+m^2)\epsilon (p_0)e^{i(p_0-m)\lambda ^2(t-t')-i{\mathbf {p}}\lambda ({\mathbf {x}}-\mathbf {x'})}}{1-e^{-\lambda ^2\beta (p_0-m+\lambda ^{-2}{\bar{\mu }})}} \end{aligned}$$

(51)

$$\begin{aligned}&\quad {\mathop {\rightarrow }\limits ^{\lambda \rightarrow \infty }}\frac{1}{(2\pi )^3}\int {d}^3{\mathbf {p}}\frac{e^{i\frac{{\mathbf {p}}^2}{2m}(t-t')-i{\mathbf {p}}({\mathbf {x}}-{\mathbf {x}}')}}{1-e^{-\beta (\frac{{\mathbf {p}}^2}{2m}+{\bar{\mu }})}}=\omega _{\beta ,{\bar{\mu }}}(\psi ^*(t,{\mathbf {x}})\psi (t',{\mathbf {x}}')) \end{aligned}$$

(52)

which is the 2-point function of the $\beta$-KMS state with chemical potential ${\bar{\mu }}$ of the nonrelativistic free scalar field $\psi$. For the second class (with condensate), we set ${\bar{\mu }}=0$ and consider a sequence of condensates

$$\begin{aligned} c_{\lambda }({\mathbf {x}})=(2m\lambda ^3)^{-\frac{1}{2}}c(\lambda ^{-1}{\mathbf {x}}) \end{aligned}$$

(53)

where c is a harmonic function. Then, the states $\omega ^{+}_{\lambda ^2\beta ,c_{\lambda }}$ converge to the state of $\psi$ with the 2-point function

$$\begin{aligned} \omega _{\beta ,c}(\psi ^*(t,{\mathbf {x}})\psi (t',{\mathbf {x}}'))=\frac{1}{(2\pi )^3}\int {d}^3{\mathbf {p}}\frac{e^{i(\frac{{\mathbf {p}}^2}{2m})(t-t')-i{\mathbf {p}}({\mathbf {x}}-{\mathbf {x}}')}}{1-e^{-\beta \frac{{\mathbf {p}}^2}{2m}}}+\overline{c({\mathbf {x}})}c({\mathbf {x}}'). \end{aligned}$$

(54)

Note that the contributions of antiparticles disappear in both cases in the limit $\lambda \rightarrow \infty$ due to the fact that the chemical potential $\mu$ tends to $+m$. Replacing $+m$ by $-m$ exchanges the role of particles and antiparticles. Note, furthermore, that the Hermitian scalar field does not have a meaningful nonrelativistic limit. Actually, one sees that the corresponding quantum states are not stable against local perturbations. This happens because local perturbations do not commute with particle number. We now pass to analyze the charge density in the nonrelativistic limit. Let $J_0d^3{\mathbf {x}}=-i:\varphi ^*{\dot{\varphi }}-\varphi {\dot{\varphi }}^*:_Hd^3{\mathbf {x}}$ denote the charge density of the complex scalar field. We scale ${\mathbf {x}}$, $\beta$ and $m-\mu$ as before and obtain

$$\begin{aligned} \lim _{\lambda \rightarrow \infty }\omega _{\lambda ^2\beta ,m-\lambda ^{-2}{\bar{\mu }}}(J_0)d^3(\lambda {\mathbf {x}})=\left( \int d^3{\mathbf {p}}\frac{1}{e^{\beta (\frac{{\mathbf {p}}^2}{2m}+{\bar{\mu }})}-1}\right) d^3{\mathbf {x}}\end{aligned}$$

(55)

for ${\bar{\mu }}>0$ and

$$\begin{aligned} \lim _{\lambda \rightarrow \infty }\omega _{\lambda ^2\beta ,c_{\lambda }}(J_0)d^3(\lambda {\mathbf {x}})=\left( \int d^3{\mathbf {p}}\frac{1}{e^{\beta \frac{{\mathbf {p}}^2}{2m}}-1}+|c({\mathbf {x}})|^2\right) d^3{\mathbf {x}}\end{aligned}$$

(56)

i.e. the charge density tends to the particle density in the nonrelativistic limit.

First Hadamard Coefficients for D

We compute in this section the first Hadamard coefficients for the operator D given in (24). We recall that in this case the Hadamard singularity has the following structure

$$\begin{aligned} H +\frac{i}{2}\Delta = \lim _{\epsilon \rightarrow 0^+}\frac{U}{\sigma _\epsilon } + V \log \left( \frac{ \sigma _\epsilon }{\xi ^2}\right) , \quad V = \sum _n V_n\sigma ^n \end{aligned}$$

(57)

similar to (2) with vanishing W, where now $U,V,V_i$ are $2\times 2$ matrices and where $\sigma _\epsilon$ is again one half of the regularized geodesic distance. In the case of Minkowski, $\sigma (x,y)=\frac{1}{2}(x-y)^\mu (x-y)_\mu$. The requirement that DH is smooth implies the following transport equations

$$\begin{aligned}&2\nabla ^\mu \sigma \nabla _\mu U + U(\square \sigma -4) -2i\mu \sigma _2 \;U \partial _0 \sigma = 0\\&- 2\nabla ^\mu \sigma \nabla _\mu V_0 - V_0(\square \sigma -2) +2i\mu \sigma _2 \;V_0\partial _0 \sigma +DU= 0\\&DV=0. \end{aligned}$$

The first two equations give

$$\begin{aligned} 2\nabla ^\mu \sigma \nabla _\mu (U^{-1}V_0) = - 2 U^{-1}V_0 - U^{-1}DU, \end{aligned}$$

considering integrals along the geodesic $\gamma$ joining x, y and indicating by r its affine parameter ($\gamma (0)=x$, $\gamma (1)=y$), that equation gives

$$\begin{aligned} 2r\frac{\mathrm{d}}{\mathrm{d} r} (U^{-1}V_0) + 2U^{-1}V_0&= - U^{-1}DU \\ 2\frac{\mathrm{d}}{\mathrm{d} r} (r U^{-1}V_0)&= - U^{-1}DU. \end{aligned}$$

Integrating along $\gamma$, as discussed e.g. in [33], we get

$$\begin{aligned} V_0(x,y) = -\frac{1}{2} U(x,y) \int _0^1 dr (U^{-1} DU). \end{aligned}$$

(58)

The first transport equation can be solved once the initial condition $U(x,x)={\mathbb {I}}$ is fixed (as we obtained in the previously computed vacuum state). To find its solution, due to translation invariance, it is enough to study

$$\begin{aligned} 2x^\mu \frac{\partial }{ \partial x^\mu } U(0,x) +2 i\mu x^0 {{\varvec{\sigma }}}_2 U(0,x) = 0 \end{aligned}$$

a solution which satisfies the desired initial condition is

$$\begin{aligned} U(0,x)= \cos (\mu x^0){\mathbb {I}} - i{{\varvec{\sigma }}}_2 \sin (\mu x^0) = \begin{pmatrix} \cos (\mu x^0) &{}\quad -\sin (\mu x^0) \\ \sin (\mu x^0) &{}\quad \cos (\mu x^0); \end{pmatrix} \end{aligned}$$

(59)

hence, since $U^{-1}= \cos (\mu x^0){\mathbb {I}} + i{{\varvec{\sigma }}}_2 \sin (\mu x^0)$

$$\begin{aligned} U^{-1}DU&= U^{-1}\left( (-\square +M^2){\mathbb {I}} + \delta M^2 {{\varvec{\sigma }}}_3 + 2i\mu {{\varvec{\sigma }}}_2 \partial _0 \right) U\\&=(M^2+\mu ^2) {\mathbb {I}} +\delta M^2 \left( \cos (2\mu x^0){{\varvec{\sigma }}}_3 -\sin (2\mu x^0) {{\varvec{\sigma }}}_1\right) \end{aligned}$$

recalling (58) we have

$$\begin{aligned} V_0(0,x)&= -\frac{1}{2} U(0,x) \int _0^1 dr\; U^{-1} DU \nonumber \\&= -\frac{1}{2} U(0,x) \left( (M^2+\mu ^2) {\mathbb {I}} +\delta M^2 \left( \frac{\sin (2\mu x^0)}{2\mu x^0}{{\varvec{\sigma }}}_3 +\frac{\cos (2\mu x^0)-1}{2\mu x^0} {{\varvec{\sigma }}}_1\right) \right) \end{aligned}$$

(60)

We can now expand V as $\sum _{n\ge 0} V_n \sigma ^n$. The equation $DV=0$ and the knowledge of $V_0$ permit to compute $V_n$ recursively. We are, in particular, interested in $[V_1](x) = V(x,x)$ because this coefficient is proportional to the trace anomaly of the stress tensor of the linearized theory [42, 55, 70] and, in particular, enters in the expressions

$$\begin{aligned} \psi _i D\psi _j, \quad \partial _a \psi _i D \psi _j. \end{aligned}$$

We observe that

$$\begin{aligned} {[}V_1] = \frac{1}{4} [DV_0]; \end{aligned}$$

furthermore, we can expand $V_0=UX$ for some X; hence, we have

$$\begin{aligned} {[}V_1] = \frac{1}{4}\left( [DU][X]-[U][\square X]-2 [\nabla _\mu U][\nabla ^\mu X]+2i\mu \sigma _2[\partial _0 X] \right) \end{aligned}$$

where

$$\begin{aligned} {[}U]={\mathbb {I}}\quad [DU] = (M^2+\mu ^2) {\mathbb {I}} + \delta M^2 {{\varvec{\sigma }}}_3 [X] = -\frac{1}{2}[DU] \\ {[}\partial _0 X] = -\frac{1}{2} \mu \delta M^2 {{\varvec{\sigma }}}_1 [\partial _0 U] = -i \mu {{\varvec{\sigma }}}_2 \quad [\square X] = -\frac{2}{3} \mu ^2 {{\varvec{\sigma }}}_3 \delta M^2. \end{aligned}$$

Summarizing this analysis

$$\begin{aligned} {[}V_1]=\frac{1}{4}[DV_0] = -\frac{1}{8} ((M^2+\mu ^2)^2+\delta M^4) {\mathbb {I}} -\frac{1}{4} \left( M^2 + \frac{\mu ^2}{3}\right) \delta M^2{{\varvec{\sigma }}}_3, \end{aligned}$$

We observe that $[V_1]$ is diagonal and constant; hence, we see that this anomaly is not visible in the conservation of the charge $J^\mu$ given in (45).

Some Technical Propositions

We report here a technical proposition, similar to Proposition 9 in [30], and adapted to the case studied here.

Proposition D.1

Consider $A_0,\ldots , A_n$ in ${\mathcal {A}}_I$, and for $k\in {\mathbb {N}}$ the following compactly supported distribution

$$\begin{aligned}&\Psi (z_1,\ldots , z_k,y_1,\ldots , y_k) \\&\quad \doteq \left. \left( \prod _{l=1}^k \frac{\delta }{\delta \psi ^{s(l)}(z_l)}\otimes \frac{\delta }{\delta \psi ^{r(l)}(y_l)} \right) (A_0 \otimes A_1\otimes \cdots A_n)\right| _{(\psi ^0,\ldots ,\psi ^n) = 0} \end{aligned}$$

where s and r maps $\{1,\ldots , k\}$ to $\{0,\ldots , n\}$ with the condition $s(l)<r(l)$. The function

$$\begin{aligned} \varphi :(p_1,\ldots , p_k)\mapsto {\hat{\Psi }}(-p_1,\ldots , -p_k,p_1,\ldots , p_k) \end{aligned}$$

given in terms of the Fourier transform ${\hat{\Psi }}$ of $\Psi$ is of rapid decrease if $P\in V_+^{k} \cup V_-^{k}$ where $V_\pm$ denotes the forward/backward light cone in the cotangent space.

Proof

Since $A_i\in {\mathcal {A}}_I$ is a microcausal functional, we have that

$$\begin{aligned} WF(\Psi (Y,Z)) \cap \left( \bigcap _{s(l)=i } W^+_l \right) \cap \left( \bigcap _{r(l)=i } W^+_{k+l}\right) = \emptyset , \quad \forall i, \end{aligned}$$

and the same holds with $W^-$ at the place of $W^+$, where

$$\begin{aligned} W^\pm _j \doteq (T^*{\mathbb {M}})^{\otimes j-1} \otimes {\overline{V}}^\pm \otimes (T^*{\mathbb {M}})^{\otimes 2k-j-1} \end{aligned}$$

for $j\in (0,\ldots , 2 k)$. Thanks to this property, we can prove that if every component of $P=(p_1,\ldots , p_k)$ is a future pointing causal vector and if $p_l = 0$ for all l such that ${r(l)=i}$ for some i, $\varphi (P)$ can be of nonrapid decrease only if

$$\begin{aligned} \sum _{{s(l)=i}} p_l = 0. \end{aligned}$$

With this observation, we can prove by induction on i that $p_l=0$ for all l; actually, if $i=0$, there are no l such that $r(l)=0$, and thus, the previous condition implies that

$$\begin{aligned} \sum _{s(l)=0} p_l = 0, \end{aligned}$$

and hence, since for all l $p_l\in V_+$, for every l such that $s(l)=0$ $p_l=0$. Furthermore, if we have already proved that $p_l=0$ for every $s(l)=j< i$, we can prove it also for $s(l)=i$. Actually, in that case we already know that for every l such that $r(l)=i$, $s(l)<i$ and thus $p_l=0$. Hence, the direction P we are analyzing can be of nonrapid decrease only if

$$\begin{aligned} \sum _{s(l)=i} p_l = 0 \end{aligned}$$

which implies that for every l such that $s(l)=i$, $p_l=0$. We have thus proved that all these directions are of rapid decrease because the direction P can be of nonrapid decrease only if $P=\{0\}$ and the zero section does not intersect $WF(\Psi )$. $\square$

In the following proposition, we prove the exponential decay of the two-point function of the KMS states for the linearized theory studied in Sect. 3.3. This proposition is used in the proof of the clustering estimate necessary for the adiabatic limit.

Proposition D.2

Let $f\in {\mathcal {E}}'(M)$, with $\text {supp} f \subset C_R$ where $C_R$ is a sphere of radius R. Consider

$$\begin{aligned} I_\sigma (u,{\mathbf {x}}) = \frac{1}{(2\pi )^3} \int d^3{\mathbf {p}} \frac{ e^{i{\mathbf {p}}{\mathbf {x}}}e^{ - \omega _\sigma u}}{(\omega _+^2-\omega _-^2)\omega _{\sigma }} \hat{{\overline{D}}}(\pm \omega _\sigma ,{\mathbf {p}}) {\hat{f}}(\omega _\sigma ,{\mathbf {p}}), \quad \sigma \in \{+,-\} \end{aligned}$$

where $\hat{{\overline{D}}}$ is as in (26), see also (24), and $\omega _\pm$ as in (25). Assume $M_->0$, it holds that

$$\begin{aligned} |I_\sigma (u,{\mathbf {x}})| \le c e^{-M_-r}, \quad r=\sqrt{|{\mathbf {x}}|^2+u^2}, \quad r>>R, \quad u>0 \end{aligned}$$

Proof

We observe that

$$\begin{aligned} -\frac{\hat{{\overline{D}}}(\pm \omega _\sigma ,{\mathbf {p}})}{(\omega _+^2 -\omega _-^2)\omega _{\sigma }}&= \left( \frac{\sigma 1}{ 2 \omega _{\sigma }} + \frac{2\mu ^2}{(\omega _+^2-\omega _-^2)\omega _{\sigma }} \right) {\mathbb {I}} + \frac{\delta M^2}{(\omega _+^2-\omega _-^2)\omega _\sigma } \sigma _3 \pm \frac{2\mu }{\omega _+^2-\omega _-^2} \sigma _2 \end{aligned}$$

Hence, we analyze separately the following functions

$$\begin{aligned} I^1_\sigma (u,{\mathbf {x}})&= \frac{1}{(2\pi )^3} \int d^3{\mathbf {p}} \frac{1}{2\omega _{\sigma }} e^{i{\mathbf {p}}{\mathbf {x}}}e^{ - \omega _\sigma u} {\hat{f}}(\omega _\sigma ,{\mathbf {p}}), \quad \sigma \in \{+,-\} \\ I^2_\sigma (u,{\mathbf {x}})&= \frac{1}{(2\pi )^3} \int d^3{\mathbf {p}} \frac{1}{2(\omega _+^2-\omega _-^2)\omega _\sigma } e^{i{\mathbf {p}}{\mathbf {x}}}e^{ - \omega _\sigma u} {\hat{f}}(\omega _\sigma ,{\mathbf {p}}), \quad \sigma \in \{+,-\} \\ I^3_\sigma (u,{\mathbf {x}})&= \frac{1}{(2\pi )^3} \int d^3{\mathbf {p}} \frac{1}{\omega _+^2-\omega _-^2} e^{i{\mathbf {p}}{\mathbf {x}}}e^{ - \omega _\sigma u} {\hat{f}}(\omega _\sigma ,{\mathbf {p}}), \quad \sigma \in \{+,-\}; \end{aligned}$$

I is then a linear combination of $I^i$ with constant coefficients. We prove with some details the decay of $I^2_\sigma$ for $\sigma \in \{+,-\}$ analogous results holds also for the other components.

Without loosing generality, let us assume that ${\mathbf {x}} = r{\mathbf {n}}$ with $n=(1,0,0)$. We shall discuss the decay for large r. Notice that

$$\begin{aligned} I^2_\sigma (u,r)&= \frac{1}{(2\pi )^3} \int d^3{\mathbf {p}} \frac{1}{(\omega _+^2-\omega _-^2)2\omega _\sigma } e^{i p_1 r}e^{ - \omega _\sigma u} {\hat{f}}(\omega _\sigma ,{\mathbf {p}}) \end{aligned}$$

We evaluate the integral in $p_1$ with the help of complex analysis considering $p_1=z=x+iy$ with $x,y\in {\mathbb {R}}$ a complex variable. We shall, furthermore, consider a particular contour of integration in the upper half plane. The integral over one branch of the contour we shall chose corresponds to $I^2_\sigma$; furthermore, the contour will be extended to infinity and chosen in such a way to avoid the poles in $1/(\omega _+^2-\omega _-^2)2\omega _\sigma$ and the brunch cuts present in $\omega _\sigma$.

Actually, we observe that since f is a compactly supported distribution, and its support is contained in the disc centered in 0 and of radius R, by the Paley Wiener theorem, its Fourier transform is an entire analytic function. Furthermore, it grows at most polynomially in every real direction and exponentially in complex directions. Hence, it exist two constants $C>0$, $C'>0$ and a $N>0$ such that

$$\begin{aligned} |{\hat{f}}(p_0,{\mathbf {p}})|\le & {} C e^{R {|\text {Im}(p_0)|+|\text {Im}({\mathbf {p}})|}} (1+|\text {Re}(p_0)|+|\text {Re}({\mathbf {p}})|)^N \nonumber \\\le & {} C' e^{R \sqrt{|p_0|^2+|{\mathbf {p}}|^2}}; \end{aligned}$$

(61)

hence, if it is composed with $\omega _\sigma$ and if it is seen as a function of $p_1$, it is analytic everywhere up the branch cuts which are present in the principal squares of

$$\begin{aligned} \omega _\pm&= \sqrt{w^2 +2\mu ^2 \pm \sqrt{4\mu ^4 + 4\mu ^2w^2 +(\delta M^2)^2}}. \end{aligned}$$

(62)

This implies that the integrand of $I^2_\sigma$, seen as a function of $p_1$, is analytic everywhere up to the poles of $1/(\omega _+^2-\omega _-^2)2\omega _\sigma$ and the branch cuts mentioned above. To describe the contour of integration, we analyze the location of the branch cuts and the poles.

We study the form of the branch cuts of $\omega _-(z,{\mathbf {p}}_\perp )$. In view of the definition of

$$\begin{aligned} \omega _-^2 = w^2 + 2\mu ^2 - 2\mu \sqrt{w^2 +\mu ^2 + \frac{\delta M^4}{4\mu ^2}} \end{aligned}$$

and recalling that $w^2 = (z^2)+|{\mathbf {p}}_\perp |^2 + M^2$ where $z=x+iy$ is the complex variables which replaces $p_1$ and ${\mathbf {p}}_\perp =(p_1,p_2)$, we have that there is a branch cut where $Z(z)\doteq w^2 +\mu ^2 + \frac{\delta M^4}{4\mu ^2}\le 0$ and where $W(z)\doteq \omega _-^2\le 0$.

The first condition, in the upper half complex plane, is met if

$$\begin{aligned} x=0, \quad y\ge \sqrt{|{\mathbf {p}}_\perp |^2 + M^2+\mu ^2 + \frac{\delta M^4}{4\mu ^2}}. \end{aligned}$$

Writing $W(z)=c+id - 2\mu \sqrt{a+ib}$, with $a,b,c,d\in {\mathbb {R}}$, the second condition is met if

$$\begin{aligned} {\left\{ \begin{array}{ll} d -2\mu \text {Im}\sqrt{Z}= d - \sqrt{2}\mu \frac{b}{\sqrt{\sqrt{a^2+b^2}+a}} =0 \\ c -2\mu \text {Re}\sqrt{Z} = c -\sqrt{2}\mu \sqrt{\sqrt{a^2+b^2}+a} \le 0 \end{array}\right. } \end{aligned}$$

(63)

In our case,

$$\begin{aligned} a&=x^2-y^2 + |{\mathbf {p}}_\perp |^2 + M^2+\mu ^2 + \left( \frac{\delta M^2}{2\mu }\right) ^2 \nonumber \\ b&= d = 2 y x \nonumber \\ c&= x^2 - y^2+ |{\mathbf {p}}_\perp |^2+ M^2+2\mu ^2 ; \end{aligned}$$

(64)

hence, since $M_->0$, $M^2 > \delta M^2 \ge 0$, (63) has solutions if

$$\begin{aligned} x=0, \quad |{\mathbf {p}}_\perp |^2 + M^2 - \delta M^2 \le y^2 \le |{\mathbf {p}}_\perp |^2 + M^2+ \delta M^2 \end{aligned}$$

or if

$$\begin{aligned} 2\mu ^2-a = \sqrt{a^2+b^2} \end{aligned}$$

which holds if

$$\begin{aligned} y = \sqrt{\frac{\mu ^2 B^2}{\mu ^2 -x^2}-\mu ^2} \end{aligned}$$

(65)

where $B^2 = |{\mathbf {p}}_\perp |^2 + M^2+\mu ^2 + \left( \frac{\delta M^2}{2\mu }\right) ^2$. Hence (65) has a solution for $|x|\le \mu$ and it holds that $y \ge \sqrt{B^2-\mu ^2}$. Furthermore, the inequality in (63) gives

$$\begin{aligned} c-2\mu ^2 \frac{b}{d} = c-2\mu ^2 \le 0 \end{aligned}$$

which gives $y^2 \ge x^2 + |{\mathbf {p}}_\perp |^2 + M^2$, which is always true.

Summarizing, the cuts are the following curves in upper half complex plane

$$\begin{aligned}&\gamma _1={\left\{ \begin{array}{ll} x=0\\ y\ge \sqrt{|{\mathbf {p}}_\perp |^2 + M^2+\mu ^2 + \frac{\delta M^4}{4\mu ^2}} \end{array}\right. } \\&\gamma _2={\left\{ \begin{array}{ll} x=0\\ \sqrt{ |{\mathbf {p}}_\perp |^2 + M^2 - \delta M^2} \le y \le \sqrt{|{\mathbf {p}}_\perp |^2 + M^2+ \delta M^2 } \end{array}\right. } \end{aligned}$$

and

$$\begin{aligned} \gamma _3={\left\{ \begin{array}{ll} |x|\le \mu \\ y = \sqrt{\frac{\mu ^2 \left( |{\mathbf {p}}_\perp |^2 + M^2+\mu ^2 + \left( \frac{\delta M^2}{2\mu }\right) ^2\right) }{\mu ^2 -x^2}-\mu ^2} \ge \sqrt{|{\mathbf {p}}_\perp |^2 + M^2 + \left( \frac{\delta M^2}{2\mu }\right) ^2}. \end{array}\right. } \end{aligned}$$

As also displayed in Fig. 1, we notice that $\gamma _1$ is contained inside of $\gamma _3$. Furthermore, $\gamma _2$ is not entirely contained in $\gamma _3$. On $\gamma _2$ and $\gamma _3$, $\omega _-$ is purely imaginary.

We observe that $\omega _+$ has only a branch cut in $\gamma _1$, furthermore, since $\omega _+^2-\omega _-^2= 4\mu \sqrt{Z}$, its branch cut is $\gamma _1$. The points where $\omega _-$ vanishes correspond to the extremal points of $\gamma _2$, those where $\omega _+$ vanishes are the extremal points of $\gamma _1$, and those where $(\omega _+^2-\omega _-^2)$ is zero the extremal points of $\gamma _1$.

Fig. 1

Tick lines correspond to the union of the brunch cuts of $\omega _-$ $\{\gamma _i\}$ while the thin line is the contour of integration $\{\xi _i\}$

The $p_1$ integration can evaluated by choosing a contour in the upper half plane. The contour can the taken in such a way that both the poles and the branch cuts lie outside it. Furthermore, the nontrivial part of the contour is chosen to be at imaginary part larger than $\sqrt{|{\mathbf {p}}_\perp |^2+M^2-\delta M^2}$. The contour $\xi$ is formed by the following paths $\xi _1,\xi _2, \xi _3, \xi _4,\xi _5$ where $\xi _1$ is the real line, $\xi _2$ is part of the semicircle centered in the origin of the complex plane and with radius B tending to infinity. Furthermore, the points with $|\text {Re}z|<2\mu$ are removed from the circle.

$$\begin{aligned} \xi _1&= \{ x \;|\; -B<x<B\}, \\ \xi _2&= \{x+iy \;|\; x^2+y^2=B^2, y\ge 0 , |x| > \mu \}, \\ \xi _3&= \{2\mu +iy \;|\; \sqrt{|{\mathbf {p}}_\perp |^2+M^2 -\delta M^2}-\epsilon<y<B\}, \\ \xi _4&= \{x +i\sqrt{|{\mathbf {p}}_\perp |^2+M^2 -\delta M^2} -i\epsilon \;|\; -2\mu<x<2\mu \}, \\ \xi _5&= \{-2\mu +iy \;|\; \sqrt{|{\mathbf {p}}_\perp |^2+M^2 -\delta M^2}-\epsilon<y<B\}. \end{aligned}$$

Notice that the poles and the branch cuts lie outside of this contour. Furthermore, in the limit $B \rightarrow 0$ the integral done on $\xi _2$ vanishes for every value of ${\mathbf {p}}_\perp$. Hence, in view of the residue theorem, and in the limit $B\rightarrow \infty$, the integral over $\xi _1$, which is nothing but $I^2_{\sigma }$, equals the sum of the integrals over $\xi _U=\xi _3\cup \xi _4\cup \xi _5$ with $B\rightarrow \infty$.

To prove that the limit $B\rightarrow \infty$ of the integral over $\xi _2$ vanishes, we observe that

$$\begin{aligned} \omega _\sigma ^2 = w^2 + 2\mu ^2 +\sigma 2\mu \sqrt{w^2 +\mu ^2 + \frac{\delta M^4}{4\mu ^2}} = \left( \sqrt{w^2 +\mu ^2 + \frac{\delta M^4}{4\mu ^2}} + \sigma \mu \right) ^2 -\frac{\delta M^4}{4\mu ^2}; \end{aligned}$$

hence

$$\begin{aligned} |\omega _\sigma |&\le \left| \sqrt{w^2 +\mu ^2 + \frac{\delta M^4}{4\mu ^2}} + \sigma \mu \right| +\frac{\delta M^2}{2\mu }\\&\le \sqrt{ \left| w^2 +\mu ^2 + \frac{\delta M^4}{4\mu ^2}\right| } + \mu +\frac{\delta M^2}{2\mu }\\&\le \sqrt{ x^2 + y^2 + |{\mathbf {p}}_\perp |^2 + M^2+ \mu ^2 + \frac{\delta M^4}{4\mu ^2}} + \mu +\frac{\delta M^2}{2\mu }\\&\le \sqrt{x^2 + y^2+ |{\mathbf {p}}_\perp |^2 + M^2-\delta M^2 } + 2\mu +\frac{\delta M^2}{\mu } +\delta M \end{aligned}$$

Furthermore, for large values of |z| if we stay outside of the region where the branch cuts are located, we have that

$$\begin{aligned} \frac{1}{|\omega _\sigma |}\le \frac{C}{|z|}, \quad \frac{1}{|\omega _+^2-\omega _-^2|}\le \frac{C}{|z|}. \end{aligned}$$

Hence, if r is sufficiently large, the exponential growth in the estimate of ${\hat{f}}(\omega ,p)$ is controlled by $|e^{ip_1r}|\le |e^{-ry }|$; furthermore,

$$\begin{aligned} \frac{1}{|\omega _\sigma ||\omega _+^2-\omega _-^2|} \end{aligned}$$

is bounded toward $\xi _2$; hence, in the limit $B\rightarrow \infty$, the integral over $\xi _2$ vanishes.

It remains to analyze the contribution of $\xi _U$. We need the following estimates. In particular, in the region contained in $\cup _i\xi _i$, the real part of any square root is positive; we have that

$$\begin{aligned} |\omega _+^2| \ge |\text {Im} \omega _+^2 | = 2|xy| \left| 1 + \frac{\sqrt{2}\mu }{\sqrt{\sqrt{a^2+b^2}+a}} \right| \ge 2|xy|, \quad z\in \xi _U \end{aligned}$$

where a and b are given in (64). For $z\in \xi _4$, since there $-y^2 + |{\mathbf {p}}_\perp |^2 + M^2>0$, it holds that

$$\begin{aligned} |\omega _+^2|\ge & {} \left| \text {Re} \left( w^2+2\mu ^2+2\mu \sqrt{w^2 +\mu ^2 +\frac{\delta M^4}{4\mu ^2}} \right) \right| \\\ge & {} 2\mu \sqrt{\left| \text {Re} \left( w^2 +\mu ^2 +\frac{\delta M^4}{4\mu ^2}\right) \right| } \ge 2\mu , \quad z\in \xi _4. \end{aligned}$$

Furthermore,

$$\begin{aligned} |\omega _+^2-\omega _-^2| \ge \sqrt{|\text {Im} (\omega _+^2-\omega _-^2)^2 |} = \sqrt{2|xy|}, \quad z\in \xi _U \end{aligned}$$

and on $\xi _4$; since there $-y^2 + |{\mathbf {p}}_\perp |^2 + M^2+2\mu ^2>0$, it holds that

$$\begin{aligned}&|\omega _+^2-\omega _-^2|^2 \ge \left| \text {Re} \left( w^2+\mu ^2+ \frac{\delta M^4}{4\mu ^2}\right) \right| = \left| x^2-y^2 + |{\mathbf {p}}_\perp |^2 + M^2+2\mu ^2\right| \ge \mu ^2, \\&\quad z\in \xi _4 \end{aligned}$$

On $\xi _4$, we have that a given in (64) is such that $a>\frac{1}{4\mu ^2}(\delta M^2+2\mu ^2)^2$; hence,

$$\begin{aligned}&|\omega _-^2|\ge |\text {Im}(\omega _-^2)| \ge 2|xy| \left| 1 - \frac{\sqrt{2}\mu }{\sqrt{\sqrt{a^2+b^2}+a}} \right| \ge 2|xy| \left| 1 - \frac{\mu }{\sqrt{a}} \right| = 2|xy| \frac{\delta M^2}{ \delta M^2 + 2\mu ^2} , \\&\quad z\in \xi _4; \end{aligned}$$

furthermore, on $\xi _3\cup \xi _5$

$$\begin{aligned} \sqrt{\sqrt{a^2+b^2}+a} \ge \sqrt{2}\mu \end{aligned}$$

and it is a monotonically decreasing function of y. Its supremum on $\xi _3\cup \xi _4$ is reached at $y^2=|{\mathbf {p}}_\perp |^2+M^2-\delta M^2 -\epsilon$ and there $a\ge \mu ^2+(\mu +\frac{\delta M^2}{2\mu })^2\ge 2\mu ^2$ and $b\ge 2\mu \sqrt{\mu ^2+(\mu +\frac{\delta M^2}{2\mu })^2} \ge 2\sqrt{2} \mu ^2$; hence,

$$\begin{aligned} \sup _{\xi _3\cup \xi _4}{\sqrt{a^2+b^2}+a} \ge (\sqrt{12}+2)\mu ^2. \end{aligned}$$

Hence, for any $\frac{1}{2}< l^2 < \sqrt{\sqrt{3}+1}-1$ we can find a ${\tilde{y}}$ where $\sqrt{\sqrt{a^2+b^2}+a}=\sqrt{2}\mu (1+l^2)$ and

$$\begin{aligned} \text {Re} (-\omega _-^2)= & {} y^2-x^2 -|{\mathbf {p}}_\perp |^2-M^2-2\mu ^2 \\&+2\mu \sqrt{2}\sqrt{\sqrt{a^2+b^2}+a} \ge y^2-x^2 -|{\mathbf {p}}_\perp |^2-M^2 +2\mu ^2 l^2, \quad y\ge {\tilde{y}} \\ \text {Re} (-\omega _-^2)\ge & {} -\mu ^2 +\delta M^2 +2\mu ^2 l^2 \ge (l^2-1) \mu + \delta M^2, \quad y\ge {\tilde{y}} \\ |\text {Im} (\omega _-^2)|= & {} |2xy| \left| 1 - \frac{\sqrt{2}\mu }{\sqrt{\sqrt{a^2+b^2}+a}}\right| \ge |2xy| \frac{l^2}{l^2+1} \quad y\le {\tilde{y}} \\ |\text {Im} (\omega _-^2)|\ge & {} |4\mu \sqrt{M^2-\delta M^2 -\epsilon }| \quad y\le {\tilde{y}} \end{aligned}$$

Combining all these estimates, we have that on $\xi _3$ and $\xi _5$

$$\begin{aligned} \left| \frac{1}{(\omega _+^2-\omega _-^2)\omega _\sigma }\right| \le C. \end{aligned}$$

uniformly in ${\mathbf {p}}_\perp$ and $\epsilon$ and the same holds true for $\frac{1}{(\omega _+^2-\omega _-^2)\omega _+}$ and for $\frac{1}{(\omega _+^2-\omega _-^2)}$ on $\xi _4$ while $\frac{1}{\omega _-}$ is bounded by an $L^1$ function uniformly in $\epsilon$ and ${\mathbf {p}}_\perp$ on $\xi _4$. Furthermore,

$$\begin{aligned} \frac{1}{|\omega _-|} \le |E|, \quad z\in \xi _4. \end{aligned}$$

With this observation, we can now control the integrals over $\xi _U$. As an example, consider the contribution to $I^2_{-}(u,r)$ due to $I^2_{-,\xi _5}(u,r)$

$$\begin{aligned} |I^2_{-,\xi _5}(u,r)|&\le \frac{1}{(2\pi )^3} \int d{\mathbf {p}}_\perp \int _{\sqrt{|{\mathbf {p}}_\perp |^2+M^2-\delta M^2}-\epsilon }^\infty \mathrm{d}y \left| \frac{1}{(\omega _+^2-\omega _-^2)\omega _-}\right| \\&\quad e^{-r y} e^{R \sqrt{2}(\sqrt{|{\mathbf {p}}_\perp |^2+M^2-\delta M^2}+|y|)} \\&\le \frac{e^{r\epsilon }}{(2\pi )^3} \int d{\mathbf {p}}_\perp \\&\quad \int _{0}^\infty \mathrm{d}y |C| e^{-(r-2\sqrt{2} R ) (\sqrt{|{\mathbf {p}}_\perp |^2+M^2-\delta M^2}+|y| )} \\&\le e^{r\epsilon }e^{-(r-2\sqrt{2} R ) \sqrt{ M^2 - \delta M^2}} \frac{1}{(2\pi )^3} \int d{\mathbf {p}}_\perp \\&\quad \int _{0}^\infty \mathrm{d}y |C| e^{-(r-2\sqrt{2} R ) (\sqrt{|{\mathbf {p}}_\perp |^2+M^2-\delta M^2}-\sqrt{ M^2 - \delta M^2}+|y|)} \\&\le e^{r\epsilon }e^{-(r-2\sqrt{2} R ) \sqrt{ M^2 - \delta M^2}} \frac{1}{(2\pi )^3} \int d{\mathbf {p}}_\perp \\&\quad \int _{0}^\infty \mathrm{d}y |C| e^{-(3-2\sqrt{2})R (\sqrt{|{\mathbf {p}}_\perp |^2+M^2-\delta M^2}-\sqrt{ M^2 - \delta M^2}+|y|)} \end{aligned}$$

where we used (61) and in the last inequality we used the fact that $r>3Rc$ and the inequalities hold uniformly in u. Hence, the integral can be taken and it can be bounded by a constant uniformly in r to get that

$$\begin{aligned} |I^2_{-,\gamma _5}(u,r)|&\le C e^{-r (\sqrt{ M^2 - \delta M^2}-\epsilon )} \end{aligned}$$

and the constant does not depend on r or $\epsilon$ for $r > 3R$. Since the inequality holds for every $\epsilon$, we have that

$$\begin{aligned} |I^2_{-,\gamma _5}(u,r)|&\le C e^{-r \sqrt{ M^2 - \delta M^2}} \end{aligned}$$

The estimate of $I^2_{-,\gamma _3}$ can be done in the same way. The analysis of $I^2_{-,\gamma _4}$ can also be done in the same way substituting |C| with the bounding $L^1$ function |E|. We actually have

$$\begin{aligned} |I^2_{-,\xi _4}(u,r)|&\le \frac{1}{(2\pi )^3} \int d{\mathbf {p}}_\perp \int _{-\mu }^\mu \mathrm{d}x \left| \frac{1}{(\omega _+^2-\omega _-^2)\omega _-}\right| e^{-r \sqrt{|{\mathbf {p}}_\perp |^2+M^2-\delta M^2}-\epsilon } \\&\quad e^{R \sqrt{2}(\sqrt{|{\mathbf {p}}_\perp |^2+M^2-\delta M^2}+|x|)} \\&\le \frac{e^{r\epsilon }}{(2\pi )^3} \int d{\mathbf {p}}_\perp \int _{-\mu }^\mu \mathrm{d}x |E| e^{-(r-2\sqrt{2} R ) (\sqrt{|{\mathbf {p}}_\perp |^2+M^2-\delta M^2}+\mu )} \\&\le e^{r\epsilon }e^{-(r-2\sqrt{2} R ) \sqrt{ M^2 - \delta M^2}} C' \\&\quad \int d{\mathbf {p}}_\perp e^{-(3-2\sqrt{2})R (\sqrt{|{\mathbf {p}}_\perp |^2+M^2-\delta M^2}-\sqrt{ M^2 - \delta M^2})} \end{aligned}$$

and it can be treated in the same way as before. All the other contributions can be analyzed in a similar way. $\square$