Finite temperature spectral functions in the O(N)-model

We directly calculate spectral functions in the O(N)-model at finite temperature within the framework of the Functional Renormalization group. Special emphasis is put on a fully numerical framework involving four-dimensional regulators preserving Euclidean O(4) and Minkowski Lorentz invariance, an important prerequisite for future applications. Pion and sigma meson spectral functions are calculated for a wide range of temperatures across the phase transition illustrating the applicability of the general framework for finite temperature applications. In addition, various aspects concerning the interplay between the Euclidean and real time two-point function are discussed.


I. INTRODUCTION
The access to real time correlation functions is key to the theoretical understanding of many interesting physics phenomena, ranging from the bound state spectrum of the theory at hand over decays to the dynamical evolution of both, close-to-and far-from-equilibrium systems. However, strongly correlated systems are typically only accessible within numerical computations with either lattice or functional approaches. These approaches have been mostly used in their Wick rotated form in Euclidean space if it comes to applications to strongly correlated systems in higher dimensions such as QCD.
In the present work we report on progress within a fully numerical approach to directly compute real time correlation functions, introduced in [1,2]. This formalism is based on the functional renormalisation group (FRG) [3], a suitable non-perturbative method, for reviews see e.g. [4,5]. Such a numerical framework is required for application to strongly-correlated systems whose degrees of freedom have a non-trivial dispersion relation. Furthermore, already in Euclidean applications the use of regulators preserving Euclidean O(4) and Minkowski Lorentz invariance is crucial in situations with finite temperatures and density. We would also like to emphasise that the current real time approach is useful and applicable to finite density situations in Euclidean space: the O(4)invariant regulators and methods discussed here allow us to guarantee the silver-blaze property without violating Lorentz symmetry. Moreover, the flow is local in both momentum-and frequency space. It is this combination of properties which is crucial for a quantitative determination of the density. The above requirements are met in QCD at finite temperature and density but also in ultracold atomic systems. Both systems are key application areas of the current real time approach.
FRG approaches to real time computations close to the present one have been set-up in [6,7]. While the former one shares the space-time symmetry with the present approach, the latter one shares the locality in momentum space. The present approach combines both momentumand frequency locality as well as Lorentz invariance.
Real-time applications range from the dynamics of low dimensional systems to the description of far-fromequilibrium dynamics, for a selection of works see  Most of the implementations have in common that either they are not directly applicable to higher dimensional systems or one would have to re-derive the plethora of Euclidean results already obtained in the Euclidean version of the approach. This situation asks for an approach which allows us to either make direct use of available Euclidean results or at least be able to benchmark the real time results with the Euclidean ones already obtained.
In addition, the physical content of Euclidean truncations in terms of scattering processes, resonances, particles and decay channels can be directly accessed from the real time correlation functions and identified with the equivalent Euclidean part.
In principle it is possible to solve the hard problem of obtaining real time correlation functions from reconstruction of Euclidean data. This tasks requires a very good resolution of the momentum and frequency dependence of the Euclidean correlation functions, for related FRG-work see [29][30][31][32][33][34][35][36][37][38][39][40][41]. Subsequently applied methods based on Bayesian Reconstruction, like the Maximum Entropy Method (MEM) [42][43][44][45][46][47] are able to produce remarkable results in many cases. For spectral reconstrunctions based on Euclidean FRG data see e.g. [48][49][50][51]. Unfortunately, these methods are mathematically not guaranteed to converge towards the correct correlation function. Other methods which do not have this problem [52][53][54], suffer from sign problems and require unachievable numerical precision for reliable results. All of these obstacles can be overcome by a direct numerical calculation.
Our long term goals are the hadron spectrum, transport coefficients and other real time observables within the framework of the fQCD collaboration [55], which aims at a quantitative first-principle description of QCD from the FRG. This stresses once more the importance of a fully numerical approach due to the sophisticated technical level already required at the Euclidean level to obtain quantitatively competitive results [36,[38][39][40]. This work is a necessary piece in order to close the gap between the calculation of Euclidean correlation functions and dynamical observables such as transport coefficients. arXiv:1711.07444v1 [hep-th] 20 Nov 2017 FIG. 1: Integration path in the complex q 0 plane at vanishing temperature These quantities can be obtained conveniently from the spectral functions as demonstrated with the shear viscosity [49,50] from reconstructed spectral functions.

II. O(N ) MODEL SPECTRAL FUNCTIONS AT FINITE TEMPERATURE
The calculation of real time observables such as spectral functions starting from a Euclidean formalism represents a difficult problem as it requires the analytic continuation from Euclidean to Minkowski signature. In this work we follow the formalism put forward in [1] where the analytic continuation is carried out on the level of the equation itself. Unfortunately the straightforward evaluation of Euclidean correlation function for complex momenta does not lead to the desired retarded correlation functions required for the extraction of the spectral function. The procedure in [1] thus requires a deformation of the integration contour into the complex plane with appropriate corrections in order to recover the desired real time observable. The following discussion is restricted to the simplified case, referred to as propagator approximation in the following, where only the Minkowski-momentum dependence of propagators is taken into account. This approximation obviously misses bound-state effects in n−point vertices that are hidden in a non-trivial Minkowski-momentum dependence of the vertices. However, such effects could still be captured by a dynamical hadronisation of the corresponding channels in the respective vertices within this approximation.
Here we present a formal solution which holds in the propagator approximation irrespective of the analytic structure of the propagators under consideration. Following [1], we start by recapitulating the case of vanishing temperature. In this approximation it is sufficient to consider a single diagram topology, namely the selfenergy contribution to the two-point function involving two three-point vertices. Here we are interested in con- Integration path in the complex q 0 plane at finite temperature tributions to the retarded 2-point functions as the latter relates directly to the spectral function. Given analytical propagators as a function of a complex momentum variable, it turns out that the contribution to the retarded correlation function is obtained by integration along a particular path in the complex momentum plane. The defining property of this path is that it represents a continuous deformation of the integration path along the real axis for vanishing Minkowski external momentum. This is illustrated in Fig. 1 for the case of two propagators which only carry a complex structures on the Minkowski axis, but this argument generalizes straightforwardly to propagators with a general complex structure. In the case of only simple poles the contribution arising from integrating along the simple contour C eucl relates to the desired contribution to the retarded correlator obtained by integrating along C ret by a closed contour integration. The difference between the two results can be written as a sum of residues of certain poles in the complex plane, which is exactly the situation covered in [1]. Turning to finite temperature, it turns out that the simple evaluation of the Matsubara sum in the analytically continued expression does again not lead to the desired contribution to the self-energy of the retarded propagator. As demonstrated in [1] for a propagator with only simple poles, the difference between the two results can be written as a sum of residues corresponding to all poles of the propagator involving the external frequency weighted by a difference of thermal occupation numbers. This generalizes straightforwardly to the situation with propagators with a general complex structure: The correction term is obtained by evaluating contour integrals surrounding all complex structures of the propagator involving the external frequency weighted with a difference of occupation numbers n(i(p 0 + q 0 )) − n(iq 0 ). This is illustrated in Fig. 2. This procedure represents a formal solution of the analytical continuation problem in the propagator approximation. The numerical implementation of this procedure turns out to be very challenging as it requires to evaluate contour integrals in the close vicinity of first-and second-order poles (as soon as the regulator term GṘG(q) is introduced). In this work we restrict ourselves to the evaluation of finite temperature spectral functions on a given LPA'+Y solution. In this case the propagators show only simple poles and the correction term can be evaluated as a sum of residue contributions instead of an actual contour integration, further details can be found in Appendix C. Technically, the retarded two-point correlator is obtained from the Euclidean correlation function: which in turn can then be used to obtain the spectral function where Z(0) is the wave function renormalization at vanishing momentum. The spectral function in (2) is renor-malization group invariant, but is not normalised to one, This comes from the normalisation with Z(0). Typically, the physical spectral function is normalised to one,

A. Truncations
In order to meet the requirement of only having to take pole corrections into account when calculating spectral functions, we restrict our effective action to the following form with the mesonic field φ = (σ, π). The effective potential, V (σ), is split into a part only depending on the O(4)- invariant ρ = 1 2 (π 2 a + σ 2 ), and a part explicitly breaking the symmetry, thus allowing for the Goldstone bosons to acquire a finite mass Note that the linear breaking term drops out of all dynamical equations, and in particular the flow equations. While this fact is hidden in an expansion about the flowing minimum, it is apparent within an expansion around a fixed expansion point in the bare field σ, see [56]. In the present work we choose the latter expansion for both stability and for a self-consistent frequency dependence, see [57]. We choose the expansion point close but above the IR minimum <σ>, as argued in [56]. We take into account terms up to φ 14 , and the convergence of the results has been checked. The truncation (5) is refered to as LPA'+Y and we additionally consider a second truncation with Z = 1, Y = 0 corresponding to the usual LPA (Local Potential Approximation) scheme. Effectively the two dressing functions Z and Y are reabsorbed into a dressing for the pion, Z π , and one for the sigma meson, Z σ .
Finite temperature is introduced via the Matsubara formalism in the usual manner.
As argued in [1], we are mostly concerned with regulators depending on the Lorentz invariant momentum configurations q 2 and is discussed at large in said reference. Further details about the regulator in the current work are given in Appendix A and a comparison with a Lorentz invariance breaking regulator is discussed in Appendix B.

III. RESULTS
We first display the temperature and momentum dependence of the spectral functions. Furthermore, we extract the pole mass and compare it to the curvature mass. Details about the numerical implementation are given in Appendix C.

A. Spectral functions at vanishing external momentum
Here we show the temperature evolution of the pion and sigma spectral function in the LPA'+Y approximation. The spectral functions feature several district structures with a clear physical interpretation, see e.g. [58] for a detailed discussion of the different processes in the Quark-Meson model. In general there are two different cut structures at finite temperature and vanishing external momentum in the propagator extended to the complex plane. The unitarity cut spanning from the multiparticle decay threshold to infinity, i.e. ω ∈ [µ thresh , ∞), which is present in any interacting theory Furthermore, the Landau cut at smaller frequencies, which is purely medium dependent and gives rise to inverse scattering processes [59] with the heat bath. These scattering processes give rise to Landau damping, hence the name. Finally, delta functions represent stable particles.
Our result for the spectral functions at vanishing external momentum in the LPA and LPA'+Y truncation are depicted in Fig. 3. For the sigma meson there are two different processes available, i.e. σ * → π + π and σ * → σ + σ and no Landau cut structure. While for the pion we have π * → π + σ for the Unitarity cut and π * + π → σ for the Landau cut. At vanishing temperature we have the expected stable pion and no stable sigma particle. At finite temperature the sigma meson emerges as stable particle as O(4)-symmetry gets restored. The presence of stable particles at finite temperatures is an artefact of the current truncation.
The difference between the two truncation is most prominently seen at multiple particle decay thresholds, which involve a sigma meson demonstrating the effect of the wave-function renormalization, shown in Fig. 5c, on the curvature mass, as in both truncations the proper pole mass is not coupled back into the system. A problem that is numerically not traceable in the current formalism, c.f. discussion in Appendix C. The cut structure of propagators is best seen in the imaginary part of the    two-point function, Im Γ (2) (ω), as a finite value translates directly into a cut of the propagator, the result is shown in Fig. 4. The finite part at small frequencies, that vanishes for larger temperatures, in the pion shows again the Landau cut. For larger frequencies the unitarity cut shows clear decay thresholds, i.e. in the sigma meson the different thresholds and their degeneracy for high temperatures can be seen nicely.
Our results compare qualitatively well to the results obtained using a spatially flat regulator, c.f. the discussion in Appendix B, in the Quark-Meson model [60][61][62][63].

B. Pole, screening & curvature masses
It is interesting to compare the real time pole and screening masses with the Euclidean curvature mass, for a detailed discussion of the respective definitions in the present FRG context see e.g. [57]. The curvature mass is typically used in Euclidean computations within low energy effective field theories for QCD. There, the physical pion and sigma masses are input parameters and are identified with the respective curvature masses. However, they have to be identified with the pole mass, and hence we have to check how well such an identification works. In the case of a stable particle, the position of the pole can be directly extracted from the spectral function. For particles with a finite decay width, this is not possible any more as the pole leaves the physical sheet of the complex energy plane and is found on the second Riemann sheet.
On the other hand the pole mass m pole = 1/ξ t is the inverse temporal screening length, which can be extracted from lim t→∞ G E (t, 0) ∝ e −t/ξt .
Furthermore the screening mass is m screen = 1/ξ spat , that is the inverse spatial screening length, and the curvature mass is Note that the wave function renormalisation Z in (9) ensures the RG-invariance of the latter. However, the necessity of a choice of the momentum configuration (in (9) it is p 0 = 0, p → 0) introduces a scheme dependence, while the pole and screening masses are schemeindependent.
Note also that at finite temperature the order of limits in (9) matters, as the two limits do not commute any more where the order of limits on the LHS is referred to as static and the limit on the RHS as plasmon, for details see e.g. [64], if not stated otherwise Z(0, 0) is understood in the static limit. The difference between the two limits can be seen easily in the corresponding spectral functions, as the static limit contains the transport peak, while this structure is absent the plasmon limit, compare e.g. Fig. 3 with Fig. 7 and the see the corresponding discussion in Section. III C. The difference between the two limits can also be seen in the corresponding Euclidean correlation functions, shown in Fig. 5b, where the green line shows the propagator as a function of the external momentum at the zeroth Matsubara mode, i.e. approaching static limit. The red line shows the propagator as a function of the external Euclidean frequency, continuously defined due to the analytic continuation involved, which approaches the plasmon limit. As soon as a small external momentum is introduced, the static limit is implicitly taken and the scale is given by the external momentum, shown as vertical lines. For better visibility the trivial mismatch for larger p has been subtracted, We can now use (9) together with (11) together with the spectral representation of the Euclidean propagator to derive a relation between the curvature mass and the spectral representation with the normalised spectral function (4). In case a stable particle is present, the spectral function can be split in two positive parts as follows, ρ =ρ pole +ρ cut , The normalisation of the pole is smaller than one: Z pole = 1 − dω 2ρ cut < 1 which follows from the normalisation of the spectral function and the positivity of ρ cut . Using this split in (11) we arrive at Eq. (14) entails the information when the difference between pole and curvature masses grows large: Firstly, decreasing Z pole increases the importance of the cut-part and hence the difference between curvature and pole mass grows. Secondly, if the spectral weight of ρ cut is taken at smaller spectral values, the integrals in (14) grow and hence the difference between curvature and pole mass grows.
Translated back to Euclidean space-time, both processes lead to strong frequency and momentum dependences in the Euclidean propagator. In turn, if the pole term dominates the full spectral function, the full Euclidean propagator is well described by a propagator with a constant wave function renormalization, depicted in Fig. 5c. A similar conclusion was also drawn in [57,65], where the relation between pole and curvature mass has also been investigated.
A very good estimate for the pole mass can be obtained from a Padé approximant of the propagator around the zero crossing of the real part of Γ (2) , if the pole is sufficiently close to the Minkowski axis. This is certainly the case for the spectral functions depicted in Fig. 3. Our result for the different masses is shown together with the order parameter σ in Fig. 5d.
While the difference is negligibly small for the pion at very small temperatures, as expected since the momentum dependence of the two-point correlator, shown in Fig. 5a, is very close to unity. A small deviation can be seen around the phase transition, where the momentum dependence becomes maximal. The more interesting case is the sigma, as its spectral function has no clear mass pole for low temperatures. Across the phase transition the curvature mass can be used as an order parameter, since it is directly related to the chiral condensate, see e.g. [56]. The pole mass of the sigma mesons behaves more gentle across the phase transition in comparison to the curvature mass, but still exhibits a clear minimum at the cross over. The larger mismatch between the two masses can again be explained by the significantly stronger momentum dependence of the sigma meson at small and medium temperatures compared to the pion, c.f. Fig. 5a. In general the qualitative strength of the mass difference can already be obtained from the temperature dependence of the constant wave function renormalizations, as Z k (0, | p|) ≈ Z k (0, k).
Furthermore, one can use (7) in order to extract the pole mass from the two-point function. Combining the spectral representation (11) with (7) one obtains i.e. the Fourier transformed propagator, which reduces to a calculation of a Laplace transform. In the case of a mass pole we were able to extract the correct pole mass from (15) again, an example for the pion at a temperature 138 MeV is shown in Fig. 6. Eq. (15) allows us to calculate the contributions from different structures in the spectral function individually, i.e. we find the expected exponential decay (7) from the mass pole and empirically the contribution from the Landau cut is very well described by additionally introducing a quadratic time dependence in the exponent. Furthermore, their sum is already sufficient to describe the full time dependence at reasonable times, for extremely large times the behaviour is trivially dominated by the necessary numerical cut η min . Unfortunately we were unable to extend this definition of the pole mass to the regime without a pole mass in the spectral function, due to a combination of a lack of the functional form of the Landau cut and numerical uncertainties.

C. Spectral functions at finite external momentum
Since there is only a very small difference between LPA and LPA'+Y in our current settings, the results with finite external momentum are calculated in LPA due to reduced numerical cost, due to which we also refrained from extrapolating our results to ε → 0. The results are depicted in Fig. 7 for various temperatures and external momenta. The main differences between spectral functions at vanishing and finite external momentum is the uniform Lorentz boost of the mass pole and the unitarity cut, as well as the transport peak, at frequencies ω < | p |, arising from space-like scattering processes at finite temperature and momenta with the heat bath, a detailed discussions about the involved kinematics can be found in [58]. Furthermore the transport peak is, like the Landau cut, not Lorentz invariant as it couples directly to the heat bath.

Lorentz invariance
A non-trivial consistency check of our results at finite external momenta is obtained by looking at the Lorentz invariance at a vanishing (or small) temperature. A Lorentz invariant function must only depend on ω 2 − p 2 , i.e. for the spectral function this translates to the property that ρ( ω 2 + p 2 , p) must be independent of p. Our results for the spectral functions at finite external momenta are depicted for this momentum configuration in Fig. 8 at a small temperature. The most notably breaking of Lorentz invariance is introduced by the finite value of the small parameter ε, but we also do not expect invariance for these parts of the spectral function. Additional breaking is visible for small frequencies in the sigma spectral function and for frequencies around the mass pole in the pion spectral function due to the small temperature present, c.f. the discussion in Section. III A and Section. III C. The remaining parts of the spectral functions, i.e. the position of the pole, the thresholds and the continuum part, show perfect Lorentz invariance. This is in contrast to most previous studies of spectral functions within the FRG where a regulator of the form (B1) was used, and therefore Lorentz invariance explicitly broken [7,58], and demonstrates the strength of our approach.

IV. CONCLUSIONS
In this work we put forward a direct calculation of finite temperature spectral functions with the FRG in the O(N ) model. This direct computation is based on a O(4)/Lorentz-invariant regularisation scheme, and can be performed on a fully numerical level. i.e. including the full momentum-and frequency dependence of correlation functions. It demonstrates the applicability of the formalism put forward in [1] at finite temperatures.
The spectral functions for the pion and sigma meson are shown across the phase transition as a function of frequency and momentum. The four expected structures in the spectral functions, i.e. the mass pole, unitarity cut, Landau cut and transport peak, are present and discussed in detail.
The spectra obtained allowed us to investigate the important relation between the curvature and the pole mass. We found that the definition of the pole mass as inverse temporal length is accessible within our framework, but unsuitable if particles are unstable and a pole is unidentifiable in the spectrum. An analytic relation between the pole and curvature mass was derived in the case of a stable particle and qualitatively verified that the difference is driven by non trivial momentum dependencies.
Furthermore, we have explicitly verified the Lorentz invariance of the spectral function. Another major advantage of the employed framework is its numerical accessibility, which makes it easily usable in more complex theories. We hope to report on applications to Yang-Mills theory and QCD in the near future.
"SFB 1225 (ISOQUANT)". It is also supported in part by the Office of Nuclear Physics in the US Department of Energy's Office of Science under Contract No. DE-AC02-05CH11231.

Appendix A: Regulator & cut-off scale
In order to use a Lorentz invariant regulator R k (p 2 ), we need to deal with additional poles/branch-cuts necessarily introduced by the non-trivial analytic structure of R k (p 2 ). While the propagator 1/Γ (2) k is restricted to have poles on the Minkowski axis, a Lorentz invariant regulator with non-trivial momentum dependence that leaves all non-analyticities of the regulated propagator 1/(Γ (2) k (q) + R k (q)) on the Minkowski axis is to date not known. As mentioned in section Section. II we use the regulator introduced in [1], i.e.
for the shape function we chose the exponential one with m = 2. The mass-like term ∆m 2 r has the effect of pushing regulator poles away from the Euclidean axis, while only having a small impact on the analytic structure of the unregulated propagator. It is parametrised as a smooth theta function where we have chosen our parameters according to Tab. I. These parameters are chosen such that the initial conditions, Γ Λ , are unchanged compared to the standard case ∆m 2 r = 0. In the vacuum it is possible to show for LPA that the regulator poles do not contribute to a certain frequency range if p 0,max is sufficiently large and the available frequency range corresponds roughly to its value. For finite temperature this is in not true anymore, as the correction factor for simple poles at position z 0 is proportional to n B (z 0 + ip 0 ) − n B (z 0 ), which is only sufficiently suppressed for |z 0 |/T 1. Therefore we have restricted our frequency range to 720 MeV, where we have checked for explicit independence of the results on the parametrisation of ∆m 2 r . The large value of ∆m 2 r , which enables us to resolve a large range of frequencies, comes with the downside Parameter p0,max α β n that the cut-off scale Λ UV must be sufficiently large. In order to still fulfil the requirement of unchanged initial conditions compared to a vanishing ∆m 2 r we have chosen Λ UV = 8.28 GeV, therefore the interpretation of our results as an effective theory of low energy QCD is strictly speaking not possible. Nevertheless our results demonstrate the applicability of the method to extract real time correlation functions from the FRG via analytic continuation and the qualitative features stay unchanged compared to the usual O(N)-model with a lower cut-off. We are unable to fix all values to their physical ones, as we observe a loss of solution similar to [66], restricting us from tuning to arbitrary IR values. Moreover we also observed a smaller and smaller range of initial values that do not lead to a break down of the numerics as we increase our truncation. Therefore we chose f π /m π = 0.93 and m σ /m π = 2.09, which fulfils the requirement m σ /m π > 2, resulting in the sigma being an unstable particle. For all truncations the initial values where chosen such that the curvature masses agree.

Appendix B: Different regulators
In addition to the regulator described in Appendix A it is also insightful to compare to the Lorentz invariance breaking regulator (B1), which has been used in most FRG related works so far This regulator has the advantage that the Matsubara sum and the following analytic continuation can be performed analytically. Therefore we compare our results obtained with the Lorentz invariant regulator (A) against the results obtained with (B1) for the LPA case and vanishing external momenta in Fig. 9. This shows that the difference between the two regulators is marginal, which can be explained by the small breaking of Lorentz invariance for p = 0, c.f. the discussion in Section. III C, and the already very similar running of the Euclidean system for the two regulators [67] in this regime. For larger temperatures the difference is dominantly driven by the difference in the condensate, leading to a general mismatch between the spectral function.
Appendix C: Numerical details We followed the general workflow outlined in [38]. The flow equations are derive using DoFun [68], traced and optimized using the FormTracer [69], utilizing FORM [70], and solved using the frgsolver, a c++ framework developed and maintained by the fQCDcollaboration [55].
In order to implement the procedure outlined in Section. II, the equations are solved as Γ (2) (ω, p) = Γ where "Correction" term refers to the correction of the occupation numbers arising implicitly from the Matsubara sum. From the study of the flow equation using a 3D regulator it is well known [1,71] that the equations reduce to delta functions for the imaginary part when the limit ε → 0 is taken by means of the Sokhotski-Plemelj theorem. Due to our Lorentz invariant regulator (A) we are however not able to resolve the Matsubara sum analytically and must take the limit numerically, therefore greatly increasing the numerical cost of the calculation as the delta functions (and derivatives thereof) have to be resolved numerically for very small values of ε. In practice the limit is performed using Richardson extrapolation with several small values of ε, the independence of the result for different sets of ε's has been checked explicitly. From the symmetries of the two-point correlator it is known that the leading term in ε of the imaginary part behaves as O(ε 1 ) and the real part behaves as O(ε 2 ), the exponents in the extrapolation are chosen accordingly.
The correction term can in principle also take branch cuts into account, as described in Section. II, and therefore allows for a self consistent treatment. With the current methods, this is numerically insurmountable as it introduces a second competing limit for the branch cut integrals.