The low-temperature behavior of the quark-meson model

We revisit the phase diagram of strong-interaction matter for the two-flavor quark-meson model using the Functional Renormalization Group. In contrast to standard mean-field calculations, an unusual phase structure is encountered at low temperatures and large quark chemical potentials. In particular, we identify a regime where the pressure decreases with increasing temperature and discuss possible reasons for this unphysical behavior.


I. INTRODUCTION
As an effective low-energy description of Quantum Chromodynamics (QCD) the quark-meson (QM) model is widely used to compute the phase diagram of strong-interaction matter. The model shares the chiral symmetry-breaking pattern of QCD and is renormalizable. When imposing spatial homogeneity a first-order chiral transition at low temperature is found in the meanfield (MF) approximation, which ends in a critical point. In this point the phase transition is of second order and lies in the same universality class as the liquid-gas transition. Such a 'chiral critical point' is searched for in dedicated experimental programs such as the beam-energy scan at RHIC.
Standard mean-field calculations in the QM model lead to a first-order phase transition line as sketched on the left-hand side of Fig. 1 where the slope, dT c /dµ c , is negative. This is consistent with the Clausius-Clapeyron relation, if one assumes that the changes in quark number density, ∆n, and entropy density, ∆s, are both positive, when going from a dilute gas phase to a liquid-like phase. The phase diagram of the QM model has also been calculated using the Functional Renormalization Group (FRG) approach which is much more powerful than mean-field theory in describing phase transitions. There it was found that the 'Lee-Wick state' of (nearly) massless quarks is avoided in the first-order transition and the sigma field remains finite across the phase boundary [1]. However, at low temperatures, the curvature of the phase boundary bends the 'wrong way' in the sense of the Clausius-Clapeyron relation as sketched on the righthand side of Fig. 1. For dµ c /dT c > 0 and ∆n > 0 the change in entropy density has to be negative, ∆s < 0.
Although relative to the Stefan-Boltzmann pressure, a decreasing pressure with increasing temperature has been observed in previous FRG studies of the phase structure of the quark-meson model [1][2][3], this phenomenon, in- In mean-field calculations the first-order line exhibits the shape shown on the left-hand side while FRG calculations find the shape on the right-hand side.
cluding possible explanations, has not been explored in the literature to our knowledge so far. Before coming to that, we briefly recall some basics on Eq. (1). Unlike the usual Clausius-Clapeyron relation for the temperature dependence of the pressure along the transition line, which is proportional to the difference of the entropy per particle on either side of the phase transition, the ∆s in Eq. (1) refers to the difference of the entropy per volume, i.e. the discontinuity of the entropy density. On one hand, in an order-disorder chiral transition one expects the entropy per particle to increase, which together with ∆n > 0 then automatically implies ∆s > 0 also, and the distinction between entropy per particle and entropy per volume does not make a qualitative difference. On the other hand, when going from a gaseous to a liquid phase, this distinction becomes important: here the entropy per particle is expected to decrease discontinuously but both signs for ∆s, the entropy-density difference, are possible depending on the size of the density jump ∆n. Only with a sufficient density increase will the entropy density also increase, and hence ∆s > 0, although the entropy per particle decreases.
In that sense, one might conclude that the mean-field result reflects the behavior of a chiral transition while the FRG result is more consistent with a gas-liquid transition to self-bound quark matter, in which ∆n is simply not large enough, in the quark-meson model. Although this is in-line with the other observations mentioned above, the pathology arises because the size of ∆s < 0 is larger than the entropy density s on the gaseous side of the first-order line so that the entropy density itself in fact turns negative across this line at low temperatures.
The same general pattern when comparing mean-field with FRG results was also observed in the parity-doublet model to describe the liquid-gas transition of nuclear matter together with a chiral transition in the high density phase [4]. No thermodynamic inconsistency arises there, however, because the discontinuity in the density at the nuclear-matter transition is always large enough when the saturation density of symmetric nuclear matter is described correctly.

II. THEORETICAL SETUP
To set the stage we briefly summarize the pertinent features of the FRG approach for the thermodynamics of the QM model as well as for real-time correlation functions. The latter will be needed to obtain the spectral properties of the in-medium sigma and pion modes which will aid the physical interpretation of the results at low temperatures T and large quark chemical potentials µ. For details of the numerical implementation we refer to our previous work [5,6].
The FRG method is based on the Wilsonian idea of integrating out fluctuations associated to different momentum scales k. At the ultraviolet (UV) scale Λ, one has to specify an ansatz for the effective average action Γ k , which will be given by the classical action of the QM model. Fluctuations from lower momentum scales are then taken into account by solving the flow equation for the effective action or equivalently the grand potential Ω(T, µ) [7,8], The infrared (IR) regulator function R k acts as a momentum-dependent mass term and suppresses fluctuations below k, Γ k denotes the second functional derivative of the effective average action w.r.t. the fields, and the trace includes a summation over internal and spacetime indices as well as an integration over the loop momenta. In the IR limit k → 0 the full quantum effective action is obtained.
To make the functional equations tractable, approximations to the effective action are needed. In the present context they amount to a derivative expansion of the effective action, which turns the flow equations into a system of partial differential equations that can be solved numerically. Usually only the lowest-order terms in the derivative expansion are kept, which is referred to as the 'local potential approximation' (LPA). In many circumstances it has proven to be reasonably accurate in the calculation of critical phenomena.
In the LPA the Euclidean effective action of the QM model (N F = 2) at a given T and µ reads Here ψ denotes the up-and down quark fields. The mesonic pion-and sigma fields are combined as φ ≡ (σ, π), h is the Yukawa coupling and cσ an explicit chiral symmetry breaking term. In the LPA only the effective potential U k flows while, e.g., renormalizations of the kinetic energy, which are of higher order in the derivative expansion, are neglected. At the UV scale Λ the effective potential is taken to be chirally unbroken, i.e.
and the UV parameters are chosen such as to reproduce a physical pion mass and pion decay constant f π in the IR limit k → 0. The IR quark mass M ψ is essentially determined by the Yukawa coupling h. The flow equation for U k can be solved numerically using a grid method (see [1] for details) and the grand potential is obtained as U k→0 ≡ Ω(T, µ). From it various thermodynamic quantities such as pressure, entropy density and quark number density can be derived in the standard way: In addition to Ω(T, µ) we will also need the real-time mesonic spectral functions ρ π,σ (ω, p). These are obtained from the Euclidean 2-point functions by the method proposed in [5,6] which entails an analytic continuation from imaginary to real energies. The flow equations for the 2point functions Γ (2) k are obtained by taking two functional derivatives of the Wetterich equation, Eq. (2), w.r.t. the fields as displayed in Fig. 2.
As initial conditions for the retarded 2-point functions at the UV scale, we use where the limit → 0 is assumed implicitly and primes denote derivatives w.r.t. the field variable φ 2 . The flow equations for the 2-point functions can then be analytically continued to real energies by the following 2-step procedure. First, the periodicity of the bosonic ¶ k G k, s and fermionic occupation numbers with respect to the discrete external (mesonic) Euclidean energy p 0 = 2nπT is exploited, Then, the Euclidean energy p 0 is replaced by a continuous real frequency ω, The spectral functions are then obtained from the retarded 2-point functions as (12) For further details on the theoretical setup and the numerical implementation we refer to [5,6].

III. RESULTS
For the results presented below the parameters listed in Tab. 1 have been used, see also [5,6]. In the present case the constituent quark mass is higher so as to shift the CEP to higher temperatures and to induce a change of curvature of the critical first-order line.  The phase diagram for the chiral order parameter σ(µ, T ) is shown in Fig. 3. The first-order line, T c (µ c ), has regions with negative and positive slope, which will be of importance in the later discussion. It is also noteworthy that there appears a small triangle-shaped region beyond the first-order line where the σ condensate seems to be enhanced. This structure is also present in the chiral limit where it is enclosed between phase transitions of first and second order, as noticed already in [1]. We note that the positive slope of the first-order line at low temperatures necessitates a decrease of the entropy density as dictated by the Clausius-Clapeyron relation in Eq. (1). This in itself is not necessarily worrisome as discussed in the introduction, and there are physical systems which exhibit such a behavior [10]. However, as shown in Fig. 4, beyond the first order phase transition we find an extended region where s is strongly negative. On the critical line this region ends just below the point where the slope dT c /dµ c changes sign (i.e. where ∆s = 0), and it decreases slowly in temperature as µ increases. By varying the parameters of the QM model and using different 3d-regulator functions we convinced ourselves that the appearance of a s < 0 region is generic and that its extent is uniquely connected to the slope of the first-order line. By calculating the entropy densities of quarks, sigma mesons and pions separately, we find that the large negative values stem from the quark contribution, while the entropy density of the mesons is always positive. In fact, the entropy density of the sigma mesons is increased in the same region where the quark contribution is negative.

IV. DISCUSSION
Regions of negative entropy density in the phase diagram are thermodynamically unacceptable and possible reasons for this behavior have to be thought of. The first that comes to mind is that this is an artifact of the choice of regulator function [11] and/or of the truncation in the derivative expansion of the effective action to lowest order. This is clearly a possibility and has to be studied in the future by going to higher orders. An example for the influence of the truncation on the phase structure can be found in [12]. What is more appealing, however, is that there are physical reasons for this problem. Several come to mind and we will discuss them in turn.

A. Color superconductivity
From extensive studies of the Nambu Jona-Lasinio (NJL) model [13] it is known that the Fermi surface of two-flavor quark matter at low temperatures can undergo a pairing transition to a (color)superconducting state. As the QM model in its physical content is very similar to the NJL model this is also expected to happen in the present case.
To elucidate this point, a schematic model study is added here [14]. When replacing the one-meson exchange diagram by a four-point vertex, where the momentum of the exchanged meson compared to the mass can be neglected, one obtains the following NJL-type Lagrangian with G α = g 2 α /(2m 2 α ) for α ∈ {σ, π}. Here we assume that the dynamical quark mass of the order of m q = 300 MeV has already been dynamically generated and the pion mass is fixed to m π = 140 MeV. The coupling constant g π is given by g π = m q /f π with the pion decay constant f π = 93 MeV. An explicit value for g σ is not needed in the following. Applying a Fierz transformation to the particle-particle channel along the lines of Ref. [13], an attractive standard 2SC channel, ψ T Ciγ 5 τ 2 λ 2 ψ , corresponding to J π = 0 + , isospin singlet, color anti-triplet, is found. In this channel the coupling constant is largest and leads to Cooper pairing. If the contribution of the sigma meson is omitted and the pion mass is kept fixed to its vacuum value this yields a coupling H s ≈ 1/140 2 MeV −2 . For a typical NJL UV cutoff of the order of Λ = 600 MeV the coupling can be estimated by H s Λ 2 ≈ 16 which is considerably larger than typical NJL values of the order of H s Λ 2 ≈ 2, cf. Ref. [13]. The conclusion is that there is some evidence for significant superconducting gaps in the quark-meson model. We therefore expect the occurrence of 2SC pairing as soon as a Fermi sphere of quarks has formed, in particular in the low-temperature and high quark-chemical potential regime beyond the phase transition where the anomalous negative entropy density is observed. A more elaborate study of the two-flavor color superconductivity within the quark-meson model based on Eliashberg-type equations for the complex and frequency-dependent pairing gaps will be presented elsewhere [15].
To properly treat pairing phenomena within the FRG the fermion propagator, from the start, has to be of the Nambu-Gorkov form, including anomalous off-diagonal terms. Their fluctuations then need to be included within the FRG in addition to those of the chiral order parameter. To our knowledge such flow equations in higherdimensional field spaces have so far been used in the quark-meson model only for finite isospin density with pion condensation [16] or in the quark-meson-diquark model to describe the competition between chiral and diquark condensate in two-color QCD [17].

B. Inhomogeneous phases
A tacit assumption in the FRG studies of the phase diagram presented above is spatial homogeneity of all phases. This may be too restrictive. Indeed, from NJL and QM model studies in the mean-field approximation [18,19] it is known that crystalline phases with inhomogeneous chiral order parameters are favored at low T and large µ over homogeneous ones. A stability analysis of the homogeneous phase should exhibit this fact. Instabilities are signaled by zero-frequency poles of the pertinent retarded real-time propagators at finite spatial critical momenta p c .
As an example we have looked at instabilities in the 'pion-direction', signaled by the inverse static pion propagator This is the famous example of the onset of 'p-wave pion condensation', studied extensively in nuclear matter [20].
To identify such critical momenta in the FRG we have considered the pion 2-point function (the inverse of D π ) to search for zeros at intermediate RG scales k, i.e. for Γ (2) k,π (ω = 0, | p c |) = 0.
k,π (ω = 0, | p|), evaluated at σ 0,k 0 (see also  Thus, we have monitored the FRG flow of the pion 2point function at finite spatial momenta obtained with the techniques described in Ref. [6]. As an example, Fig. 5 shows the real part of Γ (2) k,π (0, | p|) for µ = 400 MeV and T = 15 MeV. One observes that a zero-crossing appears at an intermediate RG scale k somewhat above 330 MeV. There the flow equations would probably start to drive the system into this new phase which is however not reachable if the constraint of homogeneity is imposed from the outside. The zero crossing may again disappear in the IR, but this result is meaningless. An analogous behavior was recently also observed as an indication of the formation of inhomogeneous phases in the propagator for collective bosonic fluctuations in the FRG description of imbalanced Fermi gases at unitarity in Ref. [21].
Based on this insight, we map the region of instability as signaled by a zero crossing in the flow of the static inverse pion propagator Γ (2) k,π (ω = 0, | p|) in the (T, µ)plane in two different ways. In the first one, we solve the flow equation for the effective potential at a given (T, µ)-point and identify its minimum at σ 0,IR in the IR. Afterwards we evaluate the flow equation for the pion 2point function Γ (2) k,π (0, | p|) on this minimum. The region where zero crossings for some momentum p = p c occur and the corresponding scales k 0 at which they occur inside this region are plotted in Fig. 6. We note that the corresponding region for the sigma 2-point function is basically the same except that it is a bit smaller. It only extends to about the k = 580 MeV line in Fig. 6. Beyond this line p-wave pion condensation therefore seems to be favored over sigma condensation, which is likely due to the residual mass difference between the two in this region. The masses of sigma and pion approach one another only gradually with increasing temperature and chemical potential as chiral symmetry gets effectively restored.
The second way to visualize the instability in the (T, µ)-plane relies on a temporary stop of the potential evolution at some intermediate scale k > k IR with the identification of the corresponding potential minimum. Subsequently, the evolution of the 2-point functions evaluated for this minimum are also stopped at this intermediate scale. This procedure is repeated for various scales k towards the IR. The region of zero crossings for some momenta and the scales at which they occur in this way are recorded in Fig. 7. This procedure has some similarity with a co-moving Taylor-expansion solution wherein the potential flow is evaluated at a scale-dependent minimum. However, as expected both procedures agree in the infrared.
The regimes in Figs. 6 and 7 provide evidence of the appearance of a spatially modulated ground state such that within these regions in the phase diagram inhomogeneous phases might be favored over homogeneous ones. The scales k 0 at the onset thereby serve as rough estimates of the scale of the modulation that might arise in the ground state, i.e. the wave number of a chiral spiral or a chiral density wave for example. The two different ways to estimate such instabilities towards a spatially modulated phase give a rough feeling for the uncertainty in the translation of an intermediate RG scale k 0 , without precise physical meaning, into a scale of modulation in the ground state.
In order to map the region of such instabilities with fluctuations beyond mean field from the FRG more pre-cisely, a truncation that allows for a spatially modulated chiral order parameter will be needed in the future.

V. SUMMARY
In FRG studies of the QCD phase diagram based on the QM model a region of negative entropy density is found at small T and large µ. It is bounded by a firstorder transition line of positive curvature, a feature that is absent in the mean-field approximation and has first been found in Ref. [1]. The region of negative entropy density which is robust against choices of parameters and 3d-regulator functions clearly points to problems with 'traditional' treatments the QM model phase diagram within the FRG. Several reasons come to mind. First, the lowest-order truncation in the derivative expansion of the effective action might be insufficient and higherorder could remedy the problem. This possibility needs to be addressed in the future. A more intriguing alternative is of physical nature in that the basic physical assumptions that go into the effective action are not fulfilled and the FRG correctly identifies this through the wrong thermodynamics. Since the negative entropy is caused in the fermion sector an obvious possibility is a Cooper instability of the Fermi surface from attractive pion and sigma exchange. Indeed it could be verified that the region of negative entropy density features (color) superconducting pairing gaps. To properly account for pairing the fermion part of the effective action has to contain Nambu-Gorkov propagators and additional flow equations for the corresponding diquark channel have to be solved. As anticipated from QM model studies of spatial inhomogeneities of the chiral order parameter in the mean-field approximation there is also the possibility that the tacit assumption of spatially homogeneous phases is incorrect. This would entail FRG calculations with spatial inhomogeneities much like as in condensed matter physics. These are clearly much more challenging numerically. Comparing Figs. 4 and 7 it seems even quite plausible that a combination of both, a pairing instability and inhomogeneities might be possible in different regimes thus leading to a rather rich phase diagram already in the quark-meson model.
To conclude, we have shown that the 'standard' treatment of chiral effective models of QCD, such as the QM model has serious problems at low temperatures and large chemical potentials within the FRG. Most likely several physical effects have to be taken into account. In addition to the possibilities discussed in this note, it is obvious that at low temperatures nucleonic degrees of freedom are of utmost importance which complicate the problem even further.