Nonperturbative fluctuation effects of charged bosonic fields:
A quark-diquark model study at nonzero density
October 01, 2025
We study the renormalization group flow of the scale-dependent effective potential of a quark-diquark model with full field dependence at nonzero chemical potential. This includes a discussion of approximations in relation to complex bosonic fields and the Silver-Blaze property. The resulting flow equation for the scale-dependent effective potential can in principle be solved down to the infrared limit. For our quark-diquark model, which may serve as a low-energy model for dense strong-interaction matter, we find that a competition between the Bardeen-Cooper-Schrieffer singularity and bosonic fluctuations can trigger a first-order phase transition at low temperatures that turns into a second-order phase transition at a tricritical point as the temperature increases.
The phase diagram of the theory of the strong interaction, QCD, at high baryon density remains a subject of ongoing research as an understanding of the details of the phase structure is ultimately required to better our understanding of, e.g., the physics underlying neutron stars. First-principles studies of various aspects of this region of the phase diagram, where QCD is expected to be in a color superconducting state, are possible with functional methods (see, e.g., Refs. [1]–[8] for early ground-breaking studies and Refs. [9]–[15] for recent developments) but remain challenging. Even more, a systematic inclusion of interaction channels of the color-superconducting type may even be relevant for an accurate first-principles localization of the critical endpoint in the QCD phase diagram [16]. In light of the aforementioned challenges, the development of effective models for high-density QCD remains of paramount importance in order to gain a more profound understanding of the mechanisms underlying dense, strong-interaction matter, see Refs. [17]–[23] for reviews. Progress in this respect has been recently made by a discussion of the renormalization of such models [24], [25], see Refs. [26]–[29] for corresponding subsequent studies.
In our present work with a clear methodological focus, we discuss the spontaneous breakdown of a global \(\text{SU}(3)\) color and \(\text{U}(1)\) vector symmetry, \(\text{SU}(3)\times\text{U}(1)_\mathrm{V}\), at nonzero temperature and large (quark) chemical potential in a QCD-inspired QDM with two massless quark flavors. In full QCD, color symmetry is a gauge symmetry and, as such, cannot be broken spontaneously [30]. Nevertheless, our QDM serves as a simple effective low-energy model for high-density QCD, in much the same way that BCS theory provides an effective description of superconductivity within electromagnetism [31]. In this model, color superconductivity is then indicated by the formation of a nonzero diquark condensate, see, e.g., Refs. [17]–[23] for reviews. To investigate the dynamical formation of the diquark condensate, we apply the fRG approach [32] which allows us to study the effect of bosonic fluctuations beyond the MFA. To be more specific, by borrowing techniques from numerical fluid dynamics [33]–[36] (see also Ref. [37]), we can compute the full field-dependent effective potential from which the diquark condensate can be deduced. We provide a detailed discussion of the derivation of the corresponding flow equation, including the absence of spurious poles in the RG flow and the role of the Silver-Blaze symmetry. For a recent overview of the range of applications of the fRG approach, from statistical mechanics and quantum many-particle systems over high-energy physics to gravity, we refer the reader to Ref. [38].
The present work is organized as follows: In 2, we introduce a QDM and discuss the relevant symmetries with a focus on the Silver-Blaze symmetry as it is relevant for the construction of truncations. After that, we briefly discuss the fRG approach including a general discussion of the construction of truncations of the effective action in 3.1. In 3.2, we then start with the MFA which neglects fluctuation effects associated with the diquark fields. The construction of truncations including such effects at leading order in the derivative expansion is presented in 3.3. Our results for the phase diagram of the QDM can be found in 4. Further applications to other theories with complex scalar fields with nonvanishing chemical potentials are demonstrated in 5. Our conclusions are presented in 6.
In this section, we introduce the phenomenological model that we mostly use for an illustration of our considerations in the present work. This model is built up from quark and diquark fields. We denote the quark fields by \(\bar{\psi}\) and \(\psi\) with \(d_\gamma = 4\) Dirac, \(N_\mathrm{f}=2\) flavor and \(N_\mathrm{c}=3\) color degrees of freedom. The quarks are assumed to form diquark pairs of the two-flavor color-superconductor (2SC) type [2]–[4], [39]–[41]. The three corresponding diquark fields \(\Delta^*_c\) and \(\Delta_c\) with \(c=1,2,3\) also carry an index referring to the fact that the diquarks carry a combination of two colors. In fact, these fields represent antisymmetric states in color space which result from the combination of two color triplets. In the following only the color index of the quark and diquark fields shall be given explicitly whereas other indices are suppressed for convenience.
For our discussion it is convenient to introduce an additional (auxiliary) complex-valued field \(\boldsymbol{\mu}\) with \(\Re(\boldsymbol{\mu})\geq0\), which is constant in spacetime. When evaluated at a real number, \(\boldsymbol{\mu}\) reduces to the standard (quark) chemical potential, and can thus be viewed as its complex extension. In \(d=4\) Euclidean spacetime dimensions, our QDM is then defined by the following action: \[\begin{align} \addtocounter{equation}{1}\label{eq:action} S[\Phi] = \int_x \biggl[&\bar{\psi}_{c} \bigl( {\rm i} \cancel{\partial} - {\rm i} \boldsymbol{\mu} \gamma^0 \bigr) \psi_{c} + \bar{\nu}^2 \Delta_{c}^* \Delta_{c} \nonumber \\ &- \bar{\psi}_{a} \gamma_{5} C \tau^2 \Delta_{c}^* \frac{{\rm i}}{2} \epsilon_{a b c} \bar{\psi}_{b} \nonumber \\ &+ \psi_{a} C \gamma_{5} \tau^2 \Delta_{c} \frac{{\rm i}}{2} \epsilon_{a b c} \psi_{b} \nonumber \biggr] \; , \end{align}\tag{1}\] where summation over repeated indices is implied, \(C\) is the charge conjugation operator, \(\tau^2\) is the second Pauli matrix living in flavor space, \(\epsilon_{a b c}\) is the Levi-Civita symbol and \(\bar{\nu}\) is a model parameter to be specified below. For convenience, we have introduced the super field \(\Phi\), \[\begin{align} \addtocounter{equation}{1}\label{eq:super-field} \Phi = (\Delta_c, \Delta^*_c, \bar\psi_c, \psi_c, \boldsymbol{\mu}) \end{align}\tag{2}\] on the left-hand side of 1 . The integral on the right-hand side is defined as \[\begin{align} \addtocounter{equation}{1}\label{eq:integral} \int_{x} \equiv \int_{0}^{\beta} \,\mathrm{d}x_0 \int_{\mathbb{R}^{d-1}} \,\mathrm{d}^{d-1}\vec{x} \, , \end{align}\tag{3}\] where \(\beta = \frac{1}{T}\) is the inverse temperature and we impose periodic and antiperiodic boundary conditions for the bosonic and fermionic fields in the compactified temporal direction, respectively.
As can be seen in 1 , the field \(\boldsymbol{\mu}\) is coupled to the quark bilinear \(\bar{\psi}_{c} \gamma^0 \psi_{c}\). Assuming that this field is constant and real-valued, it can be considered a Lagrange multiplier which allows us to tune the net baryon number density. We add that the net baryon number is a conserved quantity which results from an invariance under \(\text{U}(1)_\mathrm{V}\) transformations: \[\begin{align} \addtocounter{equation}{1}\label{eq:U4014195V-transformation} \psi_c &\mapsto u \psi_c \, , \quad &\bar{\psi}_c &\mapsto \bar{\psi}_c u^{-1} \, , \\ \Delta_c &\mapsto (u^{-1})^F \Delta_c \, , \quad &\Delta_c^* &\mapsto \Delta_c^* u^F \, , \end{align}\tag{4}\] where \(u \in \text{U}(1)_\mathrm{V}\). Since diquarks are composites of two quarks, we have \(F=2\).
At first glance, it may seem unnatural to introduce a complex-valued chemical potential \(\boldsymbol{\mu}\), since we are often ultimately interested in real-valued chemical potentials.1 However, promoting the chemical potential to a complex quantity can at least be useful for understanding certain peculiar properties of quantum field theories at low temperature. One example is the invariance of the grand canonical partition function under shifts of the chemical potential, provided the latter remains below a critical value [43]–[45]. This invariance is known as the Silver-Blaze property [46]. In practice, the promotion of the real-valued chemical potential in the standard bare action to a complex variable establishes an additional symmetry, described by the following set of transformations [43]–[45]: \[\begin{align} \addtocounter{equation}{1}\label{eq:sbs-transformation} \psi_c &\mapsto {\rm e}^{-{\rm i} \alpha x_0} \psi_c \, , \quad &\bar{\psi}_c &\mapsto \bar{\psi}_c {\rm e}^{{\rm i} \alpha x_0} \, , \\ \Delta_c &\mapsto {\rm e}^{{\rm i} F \alpha x_0} \Delta_c \, , \quad &\Delta^*_c &\mapsto \Delta^*_c {\rm e}^{-{\rm i} F \alpha x_0} \, , \\ \boldsymbol{\mu} &\mapsto \boldsymbol{\mu} - {\rm i}\alpha \, , \end{align}\tag{5}\] where \(\alpha \in \mathbb{R}\) and the associated symmetry group is \((\mathbb{R},+)\). We shall refer to this symmetry as Silver-Blaze symmetry. At nonzero temperature, one additionally requires \(\require{upgreek} \alpha=2 \uppi T n\) for \(T > 0\) with \(n \in \mathbb{Z}\) to preserve the periodicity of bosonic fields and the antiperiodicity of fermionic fields along the temporal axis.
The Silver-Blaze symmetry is a local symmetry which puts a constraint on derivatives of the fields in the temporal direction. To be specific, temporal derivatives must always come with a suitably chosen term that depends on the chemical potential [43]–[45]. In terms of our auxiliary field \(\boldsymbol{\mu}\), we have \(D_{\pm} \equiv \partial_0 \pm F \boldsymbol{\mu}\) with \(F=1\) for the quarks and \(F=2\) for the diquarks. Loosely speaking, the field \(\boldsymbol{\mu}\) acts as a single-component gauge field similar to the covariant four-potential in, e.g., quantum electrodynamics. However, the transformation of our “gauge field" \(\boldsymbol{\mu}\) does not depend on \(x_0\), see 5 . This is in contrast to the covariant four-potential which also transforms locally. The aforementioned peculiarities of quantum field theories at low \(T\) can be directly explained by the Silver-Blaze symmetry in combination with the requirement of complex analyticity in \(\boldsymbol{\mu}\), see 7.
Last but not least, the action 1 has a global \(\text{SU}(N_\mathrm{c})\) color symmetry, \[\begin{align} \addtocounter{equation}{1}\label{eq:color-transformation} \psi_b & \mapsto \left( {\rm e}^{-{\rm i}\theta^a T^a}\right)_{bc} \psi_{c} \, , \quad &\bar{\psi}_{c} &\mapsto \bar{\psi}_{b} \left( {\rm e}^{{\rm i}\theta^a T^a}\right)_{bc} \, , \\ \Delta_{b} &\mapsto \left( {\rm e}^{{\rm i}\theta^a T^a}\right)_{bc} \Delta_{c} \, , \quad &\Delta_{c}^* &\mapsto \Delta_{b}^* \left( {\rm e}^{-{\rm i}\theta^a T^a}\right)_{bc} \, , \end{align}\tag{6}\] where the \(\theta^a\) are transformation angles and the \(T^a\) are the \(\text{SU}(N_\mathrm{c})\) generators in the fundamental representation.
In the formalism one introduces an infrared regulator \(R_k\) into the path integral. This regulator is parametrized by the RG scale \(k\) that regularizes modes of the path integral with momenta \(p^2 \lesssim k^2\), whereas modes with momenta \(p^2 > k^2\) are not affected. From a Legendre transformation of this scale-dependent path integral with respect to the sources for the fields, one then obtains an RG flow equation for the scale-dependent effective average action \(\Gamma_k\), the Wetterich equation [32]: \[\begin{align} \label{eq:wetterich-equation}\addtocounter{equation}{1} \partial_t \Gamma_k [ \Phi ] = \frac{1}{2}\,\mathrm{STr} \Big[ \big( \Gamma^{(2)}_k [ \Phi ] + R_k \big)^{-1} \cdot \partial_t R_k \Big] \, . \end{align}\tag{7}\] Here, \(\mathrm{STr}\) denotes the super trace which integrates over the temporal and spatial directions, sums over all internal indices and contributes an extra minus for fermionic degrees of freedom. Furthermore, \(\partial_t = - k \partial_k\) is the RG scale derivative. The solution of the Wetterich equation, the scale-dependent effective average action \(\Gamma_k\), interpolates between a given initial condition \(\Gamma_{k=\Lambda}[\Phi]\) at the UV scale \(\Lambda\) and the full quantum effective action \(\Gamma[\Phi]\) in the IR limit \(k \to 0\), from which all physical observables can be obtained.
The UV scale \(\Lambda\) should be considered as an additional model parameter at which we fix the initial condition for \(\Gamma_k\) such that we recover the action 1 . This also implies that the regularization scheme belongs to the definition of the model. In our present work, we choose standard \(3d\) spatial regulators for the quark and diquark fields (see, e.g., Refs. [47], [48]) with corresponding Litim regulator shape functions [49], [50] and set \(\Lambda =1\,\mathrm{GeV}\). Note that, while these regulators break Lorentz invariance, they do not break any of the symmetries discussed in 2, see Refs. [45], [51], [52] and, in particular, Ref. [53] for a more detailed discussion of regulator-induced symmetry breaking. Since 7 is a highly complicated functional differential equation and therefore in general not exactly solvable, it is required to consider truncations. We define a truncation in terms of a map \({\mathbb{T}}\) which is an endomorphism on the space of all actions. This map is inserted on the right-hand side of the Wetterich equation 7 before the field derivatives are performed, which leads to the following truncated Wetterich equation: \[\begin{align} \label{eq:truncated-wetterich-equation}\addtocounter{equation}{1} \partial_t \Gamma_k [ \Phi ] = \frac{1}{2}\,\mathrm{STr} \Big[ \big( ({\mathbb{T}} \circ \Gamma_k)^{(2)} [ \Phi ] + R_k \big)^{-1} \cdot \partial_t R_k \Big] \, . \end{align}\tag{8}\] When \({\mathbb{T}}\) is the identity, we recover the exact Wetterich equation 7 . We shall refer to the action \({\mathbb{T}} \circ \Gamma_k\) as the truncated scale-dependent effective average action. We would like to emphasize that the notion of the truncation map is solely introduced to stress the difference between the solution of the truncated Wetterich equation 8 , i.e., the scale-dependent effective average action \(\Gamma_k\), and the truncated scale-dependent effective average action \({\mathbb{T}} \circ \Gamma_k\).
Typically, the truncation map \({\mathbb{T}}\) is used to reduce the complexity of a given action \(\Gamma_k\) by extracting a set of scale-dependent couplings \(\{\mathcal{O}_i(k)\}\). These couplings are then employed to construct a new (much simpler) truncated action \({\mathbb{T}}\circ \Gamma_k\). In other words, the truncation map \({\mathbb{T}}\) reduces the information that enters the right-hand side of the Wetterich equation 7 . As a consequence, the truncation map \({\mathbb{T}}\) is inherently related to the definition of projection prescriptions \(\mathrm{proj}_{\mathcal{O}_i}\) for each coupling \({\mathcal{O}_i}(k)\), i.e., \(\mathcal{O}_i(k) := \mathrm{proj}_{\mathcal{O}_i}(\Gamma_k)\). In general, the projections occurring in the truncation map should fulfill the consistency condition \[\begin{align} \addtocounter{equation}{1}\label{eq:projections-prop} \mathrm{proj}_{\mathcal{O}_i}({\mathbb{T}} \circ \Gamma_k) \stackrel{!}{=} \mathrm{proj}_{\mathcal{O}_i}(\Gamma_k) \, . \end{align}\tag{9}\] Below, we shall discuss several truncations, such as MFA and various versions of LPA for which 9 is of high relevance. Note that LPA is the lowest order in the derivative expansion and goes beyond MFA as it includes bosonic fluctuation effects of the diquark fields.
With these truncations, we investigate SSB of the \(\text{SU}(N_\mathrm{c}=3)\times \text{U}(1)_{V}\) symmetry of the QDM, leading to the SSB pattern \[\begin{align} \addtocounter{equation}{1}\label{eq:SSB-pattern} \text{SU}(3) \times \text{U}(1)_\mathrm{V} \overset{\mathrm{SSB}}{\longrightarrow} \text{SU}(2) \times \text{U}(1)\, , \end{align}\tag{10}\] i.e., five generators of \(\text{SU}(3)\times \text{U}(1)_{V}\) are spontaneously broken. From a phenomenological standpoint, this amounts to study the formation of a nonzero diquark condensate associated with the simultaneous appearance of a gap in the quark excitation spectrum, see, e.g., Refs. [21], [54] for reviews. To this end, we shall compute the scale-dependent effective potential \(U_k\) defined by the projection prescription \[\begin{align} \addtocounter{equation}{1}\label{eq:scale-dependent-effective-potential} U_k(\bar{\Delta}, \mu) = \mathrm{proj}_U(\Gamma_k) = \frac{1}{V_d} \Gamma_k[\Phi=\Phi_0] \, , \end{align}\tag{11}\] where \(V_d\) is the \(d\)-dimensional spacetime volume and we evaluate the scale-dependent effective average action at the constant field configuration \[\begin{align} \addtocounter{equation}{1}\label{eq:super-field-evaluation} \Phi_0 = (&\Delta_c=\tfrac{1}{\sqrt{2}}\bar{\Delta}\,\delta_{c 3}, \Delta^*_c=\tfrac{1}{\sqrt{2}}\bar{\Delta}\,\delta_{c 3}, \\ &\quad \bar\psi_c=0, \psi_c=0, \boldsymbol{\mu} = \mu) \, . \end{align}\tag{12}\] Here, \(\mu\) and \(\bar{\Delta}\) are real-valued numbers. This choice implies \(\Delta^*_c\Delta_c = \bar{\Delta}^2/2\). We refer to the (global) minimum \(\bar{\Delta}_\mathrm{gs}\) of the effective potential \(U = U_{k \to 0}\) as the diquark condensate, representing an order parameter for our studies of SSB below.
Let us close this subsection with a remark on the role of the Silver-Blaze symmetry in the construction of truncations. The Silver-Blaze symmetry of \(\Gamma_k\) and its potential complex analyticity in \(\boldsymbol{\mu}\) at zero temperature and low to intermediate real \(\mu\) are of high phenomenological relevance, see also 2 and App. 7. Therefore, neither of these properties of \(\Gamma_k\) should be explicitly violated by the truncation. In practice, truncations for \(\Gamma_k\) are often formulated for real \(\mu\). Therefore, they do not allow for an unambiguous assessment of the Silver-Blaze symmetry or complex analyticity in \(\boldsymbol{\mu}\). However, if a generalization of such a truncation to complex \(\boldsymbol{\mu}\) exists,2 which does not explicitly destroy complex analyticity in \(\boldsymbol{\mu}\) and respects the Silver-Blaze symmetry, then the resulting \(\Gamma_k\) will itself obey the Silver-Blaze symmetry and may be analytic in some region of the complex \(\boldsymbol{\mu}\) plane, depending on the model parameters. Only with these generalized truncations can we eventually assess whether features like the Silver-Blaze property can appear in results at \(\boldsymbol{\mu}=\mu\). In the following, we will use a complex chemical potential to construct generalized truncations, which we then discuss with respect to analyticity and Silver-Blaze symmetry. For the actual derivation of the flow equation for the scale-dependent effective potential, however, we restrict ourselves to real chemical potentials \(\boldsymbol{\mu}=\mu\). Note that the dependence of the scale-dependent effective potential on \(\Re(\boldsymbol{\mu})\) can be studied without knowledge about its dependence on the imaginary part. We rush to add a word of caution here. Even nonanalytic generalizations of truncations (e.g., truncations depending on the real and imaginary part of the chemical potential separately) can yield flow equations for the scale-dependent effective potential which agree identically with corresponding analytic generalizations at \(\boldsymbol{\mu} = \mu\). Their difference only becomes apparent for complex values of \(\boldsymbol{\mu}\). However, due to the lack of analyticity, these truncations cannot be used to assess whether features like the Silver-Blaze property are realized at real \(\boldsymbol{\mu}=\mu\).
Let us start our discussion of concrete truncations at the level of the . As discussed at the end of the previous subsection, we discuss a generalization of standard MFA to complex-valued \(\boldsymbol{\mu}\) defined by the truncation map \({\mathbb{T}}^{\mathrm{MFA}}\): \[\begin{align} \addtocounter{equation}{1}\label{eq:mf-truncation} {\mathbb{T}}^{\mathrm{MFA}} \circ \Gamma_k[\Phi] &= \int_x \biggl[\bar{\psi}_{c} \bigl( {\rm i} \cancel{\partial} - {\rm i} \boldsymbol{\mu} \gamma^0 \bigr) \psi_{c} \nonumber \\ &\qquad\quad - \bar{\psi}_{b} \gamma_5 C \tau^2 \Delta_{a}^* \frac{{\rm i}}{2} \epsilon_{a b c} \bar{\psi}_{c} \nonumber \\ &\qquad\quad + \psi_{b} C \gamma_5 \tau^2 \Delta_{a} \frac{{\rm i}}{2} \epsilon_{a b c} \psi_{c} \nonumber \biggr]\,. \end{align}\tag{13}\] This truncation map yields a \(k\)-independent functional where the structure is inherited from the QDM action 1 . It thus fulfills all symmetries of our model including the Silver-Blaze symmetry and it does not explicitly violate complex analyticity in \(\boldsymbol{\mu}\) since the real and imaginary parts do not occur separately. However, this map removes the back coupling of the scale-dependent effective potential \(U_k\) on the RG flow as well as the effect of kinetic terms of the diquark fields. Still, among infinitely many other terms, a nontrivial scale-dependent effective potential as well as kinetic terms for the diquark fields are dynamically generated in the RG flow. For the flow equation of the scale-dependent effective potential we obtain3 \[\begin{align} \addtocounter{equation}{1}\label{eq:MF-medium} \partial_t U_k &= v_{d-1} N_\mathrm{f} d_\gamma \,k^d \sum_{\chi = \pm 1} \frac{k+\chi \mu}{2E^\mathrm{F}_\chi}\bigl(1-2 n_\mathrm{F}(\beta E^\mathrm{F}_\chi)\bigr) \, , \end{align}\tag{14}\] where \(\require{upgreek} v_{d-1} = \mathrm{vol}({d-1})/(2\uppi)^{d-1}\), \(\mathrm{vol}(d-1)\) is the volume of a (\(d-1\))-dimensional unit ball, \(n_{\mathrm{F}}(x) = 1/({\rm e}^x+1)\) is the Fermi-Dirac distribution function and \[\begin{align} \addtocounter{equation}{1}\label{eq:E-F} E^\mathrm{F}_\chi = \sqrt{\frac{\bar{\Delta}^2}{2} + (k+\chi\mu)^2} \, \end{align}\tag{15}\] are the quark energy functions.
The RG flow of the scale-dependent effective potential is predominantly driven by the so-called BCS singularity in the limit of vanishing temperature, see, e.g., Ref. [45] for a detailed discussion in the context of RG flows. More precisely, the curvature of the scale-dependent effective potential is divergent at \(\bar{\Delta}=0\) from \(k=\mu\) down to the IR limit associated with \(k \to 0\). This phenomenon is of high relevance for the numerical treatment of the QDM. Therefore, we shall discuss it in some detail in the following.
Starting with the zero-temperature limit of 14 , we find \[\begin{align} \addtocounter{equation}{1}\label{eq:MF-T0} \partial_t U_k &= v_{d-1} N_\mathrm{f} d_\gamma\, k^d \sum_{\chi = \pm 1} \frac{k+\chi \mu}{2E^\mathrm{F}_\chi}\,. \end{align}\tag{16}\] An integration with respect to the RG scale \(k\) from \(k^\prime\) to \(\Lambda>\mu\) then yields
\[\begin{align} \addtocounter{equation}{1}\label{eq:MF-pot} U_{k^\prime} &=\frac{1}{2}\bar{\nu}^2 \bar{\Delta}^2 + v_{d-1} N_\mathrm{f} d_\gamma \int^{\Lambda}_{k'} \frac{\,\mathrm{d}k}{k} \, k^{d} \sum_{\chi = \pm 1} \frac{k+\chi \mu}{2E^\mathrm{F}_\chi} \\ &=\frac{1}{2}\bar{\nu}^2 \bar{\Delta}^2 - \frac{v_{d-1} N_\mathrm{f} d_\gamma}{64} \biggl[ -2\mu(13 \bar{\Delta}^2 -4\mu^2-4k^2) \sum_{\chi = \pm 1} \chi \sqrt{\tfrac{\bar{\Delta}^2}{2} + (k+\chi\mu)^2} \\ &\quad\quad+ 2k( -4k^2+3\bar{\Delta}^2-4\mu^2 ) \sum_{\chi = \pm 1} \sqrt{\tfrac{\bar{\Delta}^2}{2} + (k+\chi\mu)^2} -3\bar{\Delta}^2(\bar{\Delta}^2-8\mu^2) \sum_{\chi = \pm 1} \mathrm{tanh}^{-1}\biggl\{ \frac{k+\chi \mu}{\sqrt{\tfrac{\bar{\Delta}^2}{2} + (k+\chi\mu)^2}} \biggr\}\biggr]_{k^\prime}^\Lambda\,, \end{align}\tag{17}\]
where we have used the action 1 . Expanding this result for \(k^\prime>\mu\) about \(\bar{\Delta} =0\) reveals that it is analytic at this point and yields a convergent Taylor series. In contrast, for \(k^\prime \leq \mu\) the result is nonanalytic at \(\bar{\Delta} =0\) due to logarithmic contributions originating from the inverse of the hyperbolic tangent, namely \(\tfrac{3}{8}\mu^2\bar{\Delta}^2\ln \tfrac{\bar{\Delta}^2}{\Lambda^2}\) for \(k^\prime<\mu\) and \(\tfrac{3}{16}\mu^2\bar{\Delta}^2\ln \tfrac{\bar{\Delta}^2}{\Lambda^2}\) for \(k^\prime=\mu\), see also, e.g., Ref. [13]. It directly follows that the curvature \(\partial^2_{\bar{\Delta}} U_{k^\prime}\) for \(k^\prime \leq \mu\) scales as \(\sim\mu^2 \ln \tfrac{\bar{\Delta}^2}{\Lambda^2}\) for \((\bar{\Delta}/\Lambda)^2 \to 0\), i.e., it tends to \(-\infty\) in this limit.4 Consequently, for any \(\mu\neq 0\), the ground state of our model is governed by SSB, i.e., a nonzero diquark condensate \(\bar{\Delta}_\mathrm{gs}\) is necessarily generated for \(k' \leq \mu\) and persists in the IR limit, irrespective of our choice for the finite parameter \(\bar{\nu}^2>0\) in the action 1 . Moreover, from an expansion of \(\partial_{\bar{\Delta}} U\) in small \((\mu/\Lambda)^2\) and \((\bar{\Delta}/\Lambda)^2\), we can deduce the characteristic BCS-like \(\mu\)-dependence of the condensate from the RG flow, see, e.g., Ref. [45]: \[\begin{align} \addtocounter{equation}{1}\label{eq:BCS-scaling} \bar{\Delta}_\mathrm{gs} \propto \displaystyle{\rm e}^{-c/\mu^2} \, , \end{align}\tag{18}\] where \(c>0\) is a \(\bar{\nu}^2\)-dependent constant.
To illustrate the manifestation of the BCS singularity in the RG flow, we show the derivative of the scale-dependent effective potential \(\partial_{\bar{\Delta}} U_k\) for \(\mu = 0.5\,\mathrm{GeV}\) as a function of \(\bar{\Delta}\) for different values of \(k\) in 1. For vanishing temperature, we find that \(\partial_{\bar{\Delta}} U_k\) is discontinuous at \(\bar{\Delta} = 0\) for \(k = \mu\) and consequently the curvature diverges. For \(k<\mu\), the discontinuity disappears but the curvature at \(\bar{\Delta} = 0\) still diverges.
At nonzero temperature, the BCS singularity is screened. This is demonstrated in 1 for \(T=5\,\mathrm{MeV}\). To be specific, we observe that the discontinuity in \(\partial_{\bar{\Delta}} U_k\) at \(\bar{\Delta}=0\) for \(k = \mu\) disappears and the region around the origin turns into a “region of steep descent", being a relic of the BCS singularity. Outside of this region, we find that our results for \(\partial_{\bar{\Delta}} U_k\) for \(T = 0\) and \(T=5\,\mathrm{MeV}\) agree almost identically in the IR limit. Indeed, we observe pointwise convergence of our finite-temperature results for \(\partial_{\bar{\Delta}} U_k\) to the zero-temperature result as the temperature is decreased. In particular, the position of the nontrivial minimum of the potential \(U_k\) (zero of \(\partial_{\bar{\Delta}} U_k\)) approaches smoothly the value obtained in the zero-temperature calculation.
We now discuss the construction of truncations that incorporate bosonic fluctuations of the diquark fields at lowest order in the derivative expansion, commonly referred to as LPAs. Due to the complex nature of the diquark fields and their coupling to the “gauge field" \(\boldsymbol{\mu}\), this task is nontrivial.
In LPA the change of the scale-dependent effective potential is taken into account in the truncation map. To enable a discussion of the Silver-Blaze symmetry, the potential – being explicitly dependent on the chemical potential – must be consistently generalized to the complex-valued chemical potential \(\boldsymbol{\mu}\). Moreover, the truncation map should preserve complex analyticity in \(\boldsymbol{\mu}\). A natural choice that does not break analyticity explicitly is \[\begin{align} \addtocounter{equation}{1}\label{eq:new-potential} \boldsymbol{U}_k(\bar{\Delta},\boldsymbol{\mu}) = \mathrm{proj}_{\boldsymbol{U}}(\Gamma_k) = \frac{1}{V_d} \Gamma_k[\Phi=\tilde{\Phi}_0] \, , \end{align}\tag{19}\] where the field configuration is now \[\begin{align} \addtocounter{equation}{1}\label{eq:new-super-field-evaluation} \tilde{\Phi}_0 = (&\Delta_c=\tfrac{1}{\sqrt{2}}\bar{\Delta}\,\delta_{c 3}, \Delta^*_c=\tfrac{1}{\sqrt{2}}\bar{\Delta}\,\delta_{c 3}, \\ &\quad \bar\psi_c=0, \psi_c=0, \boldsymbol{\mu}) \,. \end{align}\tag{20}\] Here, \(\boldsymbol{\mu}\) is still a complex number and is not yet evaluated at a real number \(\mu\) as done in 12 . The scale-dependent effective potential 11 and the generalized scale-dependent effective potential 19 are related by \[\begin{align} \addtocounter{equation}{1}\label{eq:relation-potentials} U_k(\bar{\Delta},\mu)=\boldsymbol{U}_k(\bar{\Delta},\boldsymbol{\mu} = \mu) \, . \end{align}\tag{21}\] Note that, in general, \(\boldsymbol{U}_k(\bar{\Delta},\boldsymbol{\mu})\) may depend on \(\Re(\boldsymbol{\mu})\) and \(\Im(\boldsymbol{\mu})\) separately, if it is nonanalytic in \(\boldsymbol{\mu}\).
With \(\boldsymbol{U}_k\) at hand, we can now define a quite general abstract truncation map for LPA: \[\begin{align} {\mathbb{T}}^{\mathrm{LPA}} \circ \Gamma_k[\Phi] &= {\mathbb{T}}^{\mathrm{MFA}} \circ \Gamma_k[\Phi] + \int_x \boldsymbol{U}_k(\sqrt{2\Delta^*_c \Delta_c}, \boldsymbol{\mu}) \\ &\quad + \int_x \vec{\nabla} \Delta^*_c \vec{\nabla} \Delta_c + (\text{t.d.k.}) \, . \addtocounter{equation}{1}\label{eq:LPA-truncation-map} \end{align}\tag{22}\] Here, \({\mathbb{T}}^{\mathrm{MFA}} \circ \Gamma_k[\Phi]\) is the MFA truncation map given in 13 and “t.d.k.” stands for “temporal diquark kinetic” terms (i.e., contributions involving temporal derivatives of the diquark fields), which are assumed to be scale-independent within LPA. Furthermore, we have promoted the field-argument \(\bar{\Delta}\) of the generalized scale-dependent effective potential \(\boldsymbol{U}_k\) to a position-dependent quantity, i.e., \(\sqrt{2\Delta^*_c(x) \Delta_c(x)}\). We will come back to the discussion of the Silver-Blaze symmetry in LPA in a moment, but first let us discuss possible choices of temporal kinetic terms and their implications.
As a consequence of our discussion in 2, it may be tempting to employ the following ansatz as a first version for the temporal diquark kinetic term: \[\begin{align} \addtocounter{equation}{1}\label{eq:kinetic-v1} (\text{t.d.k.})_{\mathrm{1}} = \int_x D_{-} \Delta^*_c D_{+} \Delta_c\,. \end{align}\tag{23}\] We shall refer to the truncation map with this temporal kinetic term as \({\mathbb{T}}^{\mathrm{LPA}_1}\). However, this choice leads to double counting in our calculations since the covariant temporal derivatives contribute a term \(\sim \boldsymbol{\mu}^2 \Delta^*_c \Delta_c\) to the truncation map which is already included in corresponding second-order terms in the generalized scale-dependent effective potential \(\boldsymbol{U}_k\). More explicitly, when evaluating 22 at \(\tilde{\Phi}_0\) with 23 as the temporal kinetic term for the diquarks, we obtain \[\begin{align} \addtocounter{equation}{1}\label{eq:double-counting} \frac{1}{V_d}({\mathbb{T}}^{\mathrm{LPA}_1} \circ \Gamma_k)[\tilde{\Phi}_0]&=\boldsymbol{U}_k(\bar{\Delta}, \boldsymbol{\mu}) - \frac{1}{2}F^2 \boldsymbol{\mu}^2 \bar{\Delta}^2\\ & \neq \frac{1}{V_d}\Gamma_k[\tilde{\Phi}_0] \stackrel{\eqref{eq:new-potential}}{=} \boldsymbol{U}_k(\bar{\Delta}, \boldsymbol{\mu}) \, . \end{align}\tag{24}\] For nonzero \(\boldsymbol{\mu}\), this is in contradiction to 9 ; the latter requires to leave the generalized scale-dependent effective potential unchanged. Below, in our discussion of the flow equations, which follow from our truncation maps, we shall see that this double-counting issue also leads to problems with respect to the phenomenological interpretation of the results.
From these considerations we deduce that a proper version of LPA must exclude potential-like contributions in the temporal diquark kinetic term. This issue may be solved by employing \[\begin{align} \addtocounter{equation}{1}\label{eq:kinetic-v2} (\text{t.d.k.})_{\mathrm{2}} = \int_x \Big(&D_{-} \Delta^*_c D_{+} \Delta_c + F^2 \boldsymbol{\mu}^2 \Delta^*_c\Delta_c \Big) \,. \end{align}\tag{25}\] In the following we shall refer to the truncation map with this temporal kinetic term as \({\mathbb{T}}^{\mathrm{LPA}_2}\). Note that, although not explicitly present in the truncation map, a term \(\sim \boldsymbol{\mu}^2 \Delta^*_c \Delta_c\) is dynamically generated in the RG flow [24] and eventually appears implicitly as part of the generalized scale-dependent effective potential \(\boldsymbol{U}_k(\bar{\Delta}, \boldsymbol{\mu})\) after performing the projection 11 . In our discussion below, we shall carefully examine the consequences of these two choices for the LPA truncation map, also from a phenomenological standpoint.
Let us now discuss the Silver-Blaze symmetry and whether it is realized or not for our two versions of LPA. At first glance, by looking at the temporal diquark kinetic terms, we note that \((\text{t.d.k.})_{\mathrm{1}}\) is invariant under the Silver-Blaze transformations 22 whereas \((\text{t.d.k.})_{\mathrm{2}}\) is not. However, this does not directly imply that \({\mathbb{T}}^{\mathrm{LPA}_1}\) is a Silver-Blaze-invariant truncation, since the contribution from the generalized scale-dependent effective potential must also be taken into account. In fact, both truncations are likely to break the Silver-Blaze symmetry. This is because the generalized scale-dependent effective potential in the truncation map 22 can develop a highly nontrivial dependence on \(\Im(\boldsymbol{\mu})\) during the RG flow, see also App. 7.5 Its transformation behavior under Silver-Blaze transformations can therefore, in general, not be compensated by the corresponding transformation of the simple temporal diquark kinetic terms \((\text{t.d.k.})_{\mathrm{1}}\) and \((\text{t.d.k.})_{\mathrm{2}}\), respectively.6 Moreover, depending on the initial condition for the generalized scale-dependent effective potential, the Silver-Blaze symmetry may already be explicitly broken at the UV scale.
In general, if we would like to construct a Silver-Blaze-symmetric version of LPA and at the same time would like to use the full \(\boldsymbol{\mu}\)-dependence of the generalized scale-dependent effective potential as shown in 22 , then this would require a highly complicated and RG-scale dependent temporal kinetic term for the diquark fields (in general involving derivatives of arbitrarily higher orders). Therefore, with respect to studies in LPA, we discard the Silver-Blaze symmetry and will use \({\mathbb{T}}^{\mathrm{LPA}_2}\) for all numerical calculations presented below.
For completeness, let us now discuss the special case where \(\boldsymbol{U}_k(\bar{\Delta},\boldsymbol{\mu})\) is an analytic function in \(\boldsymbol{\mu}\) in a neighborhood of \(\boldsymbol{\mu}=0\). Complex analyticity implies that \(\boldsymbol{U}_k(\bar{\Delta},\boldsymbol{\mu})\) can be written in a Taylor series in \(\boldsymbol{\mu}\) around \(\boldsymbol{\mu}=0\) and, in particular, we are allowed to promote \(\mu\) in \(U_k(\bar{\Delta}, \mu)\) to a complex number, i.e., we have \(\boldsymbol{U}_k(\bar{\Delta}, \boldsymbol{\mu}) \equiv U_k(\bar{\Delta}, \boldsymbol{\mu})\). Furthermore, in case of a complex analytic and Silver-Blaze-symmetric scale-dependent effective average action, any \(\boldsymbol{\mu}\)-dependence in the pure diquark sector has to be encoded in terms of the following building blocks: \[\begin{align} \addtocounter{equation}{1}\label{eq:building-block} D^n_{-} \Delta^*_c D^m_{+} \Delta_c \, , \end{align}\tag{26}\] where \(m\) and \(n\) are integers with \(m, n \geq 0\). Consequently, in a Silver-Blaze-symmetric version of the truncation map 22 , the pure diquark sector can be written as \[\begin{align} \addtocounter{equation}{1}\label{eq:SBS-improved-truncation} &\int_x \Big( \boldsymbol{U}_k(\sqrt{2\Delta^*_c \Delta_c}, \boldsymbol{\mu}) + (\text{t.d.k.}) \Big)\\ &\quad \stackrel{!}{=} \int_x \Big(U_k(\sqrt{2\Delta^*_c \Delta_c}, \mu=0) + g_k(\{ D^n_{-} \Delta^*_c D^m_{+} \Delta_c\}) \Big)\, . \end{align}\tag{27}\] Here, \(U_k(\sqrt{2 \Delta^*_c \Delta_c}, \mu=0)\) is the scale-dependent effective potential evaluated at vanishing chemical potential. The RG-scale dependent function \(g_k\) does not carry an explicit dependence on \(\boldsymbol{\mu}\). It is solely constructed from the building blocks given in 26 . Apparently, the right-hand side of 27 is invariant under the Silver-Blaze transformations in 5 . However, similar to the previous discussion without limitation to complex analyticity, we can deduce from 27 that the temporal diquark kinetic term is in general given by a highly nontrivial RG-scale dependent function. Consequently, in the analytic case, any truncation that includes only a finite number of temporal derivative terms of the diquark fields but the full \(\mu\)-dependence of the scale-dependent effective potential as in 22 will in general break the Silver-Blaze symmetry. However, if the truncation map includes only a Taylor expansion of the scale-dependent effective potential in \(\boldsymbol{\mu}\) up to a finite order, i.e., if the full \(\boldsymbol{\mu}\)-dependence is not taken into account, it is still possible to preserve the Silver-Blaze symmetry. Using 27 as starting point for the construction of such a truncation, it is reasonable to define the lowest order of the derivative expansion by the truncation map \[\begin{align} \addtocounter{equation}{1}\label{eq:LPA-SBS-truncation-map} {\mathbb{T}}^{\mathrm{LPA},\mathrm{SBS}} \circ \Gamma_k[\Phi] &= {\mathbb{T}}^{\mathrm{MFA}} \circ \Gamma_k[\Phi] \\ &\quad + \int_x U_k(\sqrt{2\Delta^*_c \Delta_c}, \mu=0) \\ &\quad + \int_x \vec{\nabla} \Delta^*_c \vec{\nabla} \Delta_c + D_{-} \Delta^*_c D_{+} \Delta_c \, , \end{align}\tag{28}\] where SBS stands for Silver-Blaze symmetric. This map includes \(\boldsymbol{\mu}\) dependence up to \(\mathcal{O}(\boldsymbol{\mu}^2)\). Note that, arbitrarily high orders of \(\boldsymbol{\mu}\) and derivatives of diquark fields are still generated in the RG flow. The scale-dependent effective average action \(\Gamma_k\) resulting from this map will be Silver-Blaze symmetric as well, provided that the regulator does not break the Silver-Blaze symmetry. However, since the complex analyticity of the scale-dependent effective average action cannot be guaranteed and the flow equation for this truncation exhibits systematic issues in the case of SSB, as we shall discuss next in 3.3.1, we do not employ it in our studies of specific models in Secs. 4 and 5. Nevertheless, this truncation may still prove useful in applications where the diquark field is only used as an auxiliary field to resolve momentum dependences of four-quark correlators in phases without (color-)superconducting condensate, e.g., at high temperature and chemical potential in case of QCD.
For all versions of LPA discussed above, the flow equation for the scale-dependent effective potential assumes the same functional form: \[\begin{align} \addtocounter{equation}{1}\label{eq:LPA-flow-equation} \partial_t U_k &= -v_{d-1} k^{d} \, 2(N_\mathrm{c} - 1) G_1(E^{\mathrm{B}}_1) \\ &\quad - v_{d-1} k^{d}\sum_{\chi=\pm1}\Big(1 + \chi\frac{4 F^2 \mu^2}{\xi^2_+ - \xi^2_-}\Big) G_2(\xi_\chi)\\ &\quad+ v_{d-1} k^d N_\mathrm{f} d_\gamma \sum_{\chi = \pm 1} \frac{k+\chi \mu}{2E^\mathrm{F}_\chi}\bigl(1-2 n_\mathrm{F}(\beta E^\mathrm{F}_\chi)\bigr) \, , \end{align}\tag{29}\] where \(E^\mathrm{F}_\chi\) is the quark energy function defined in 15 and \[\begin{align} \addtocounter{equation}{1}\label{eq:xis} \xi^2_{\pm} &= E_+^2 + 2 F^2\mu^2 \\ &\quad\pm \sqrt{(E^2_+ + 2 F^2\mu^2)^2 - E^4_+ + E_-^4} \end{align}\tag{30}\] with \(E^2_\pm = \frac{1}{2}(E_2^{\mathrm{B}})^2 \pm \frac{1}{2}(E_1^{\mathrm{B}})^2\). The diquark energy functions \(E^{\mathrm{B}}_1\) and \(E^{\mathrm{B}}_2\) differ for our different LPA versions. We shall define them below.
In the flow equation 29 two diquark propagator functions appear: \[\begin{align} \addtocounter{equation}{1}\label{eq:bosonic-propagator} G_1(E) &= \frac{k}{2\sqrt{F^2 \mu^2 + E^2}} \Big( 1 +\\ &\quad + n_{\mathrm{B}}(\beta (\sqrt{F^2 \mu^2 + E^2}-F \mu)) \\ &\quad + n_{\mathrm{B}}(\beta (\sqrt{F^2 \mu^2 + E^2}+F \mu)) \Big) \end{align}\tag{31}\] and \[\begin{align} \addtocounter{equation}{1}\label{eq:bosonic-propagator-2} G_2(E) &= \frac{k}{2 E} \Big( 1 + 2n_{\mathrm{B}}(\beta E) \Big) \, , \end{align}\tag{32}\] where \(n_\mathrm{B}(x)=1/({\rm e}^x - 1)\) is the Bose-Einstein distribution function. Note that \(0 \leq 4F^2\mu^2/(\xi^2_+ - \xi^2_-) \leq 1\) for all \(\mu\) and \((E^{\mathrm{B}}_1)^2>0, (E^{\mathrm{B}}_2)^2>0\).
In the limit of vanishing chemical potential the flow equation reduces to \[\begin{align} \partial_t U_k &= -v_{d-1} k^{d} (2N_\mathrm{c} - 1) G_2(E^{\mathrm{B}}_1) - v_{d-1} k^{d} G_2(E^{\mathrm{B}}_2)\\ &\quad+ v_{d-1} k^{d+1} N_\mathrm{f} d_\gamma \frac{1-2 n_\mathrm{F}(\beta E^\mathrm{F}_0)}{E^\mathrm{F}_0} \, ,\addtocounter{equation}{1}\label{eq:LPA-flow-equation-mu610} \end{align}\tag{33}\] where we have used \(\xi^2_{+}|_{\mu = 0} = (E^{\mathrm{B}}_2)^2\), \(\xi^2_{-}|_{\mu = 0} = (E^{\mathrm{B}}_1)^2\), and the fact that \(4F^2\mu^2/(\xi^2_+ - \xi^2_-)\) vanishes at least linearly in \(\mu\). The bosonic contribution can be split into two parts: one depending solely on \(E^{\mathrm{B}}_1\) and the other depending solely on \(E^{\mathrm{B}}_2\). The factor \(2 N_\mathrm{c} - 1 = 5\) in the \(E^{\mathrm{B}}_1\) contribution is the number of broken generators as determined by the symmetry-breaking pattern 10 . At nonzero chemical potential, see 29 , the bosonic contributions to the RG flow of the scale-dependent effective potential cannot be split into the aforementioned two parts, see, e.g., Refs. [56]–[58] for a more detailed analysis of this fact. Note that the mixing of the \(E^{\mathrm{B}}_1\) and \(E^{\mathrm{B}}_2\) dependence in the bosonic contribution at nonzero \(\mu\) renders the flow equation 29 numerically more challenging than, e.g., the LPA flow equation of the quark-meson model where the corresponding energy functions associated with mesons do not mix, even at nonzero \(\mu\), see App. 8 for details of the numerical implementation.
Provided that the squared diquark energy functions are greater than zero, \((E_1^{\mathrm{B}})^2 > 0\) and \((E_2^{\mathrm{B}})^2 > 0\), we have \(\xi^2_{\pm} > 0\) for all \(\mu\), see 30 . If one of the (squared) energy functions approaches zero, we have \(\xi_- \to 0\) while \(\xi_+\) remains positive, \(\xi_+ > 0\). In this case, the diquark propagator function \(G_2(\xi_-)\) diverges in the RG flow, see 29 . To be specific, we have \(G_2(E) \sim k/(\beta E^2)\) for \(E \to 0\). This behavior of the propagator prevents the RG flow from actually reaching the singularity of \(G_2(E)\) for \(k>0\). For the diquark energy functions, this implies that \[\begin{align} \addtocounter{equation}{1}\label{eq:constraints} (E_1^{\mathrm{B}})^2 \geq 0 \,, \quad (E_2^{\mathrm{B}})^2 \geq 0 \, . \end{align}\tag{34}\] for any \(k\geq 0\). In fact, as we shall discuss below, this “self-healing” property of the flow equation eventually ensures that the scale-dependent effective potential becomes convex in the IR limit, see Refs. [59], [60]. Note that the singular behavior becomes weaker as the temperature is decreased. In fact, the propagator function only diverges as \(G_2(E) \sim k / (2 E)\) at \(T=0\).
The situation is different for the function \(G_1\). While it also diverges as \(G_1(E) \sim k/(\beta E^2)\) in the limit \(E\to 0\) for \(T>0\), it converges pointwise to \(1/(2 \sqrt{F^2 \mu^2 + E^2})\) for all \(E > 0\) in the zero-temperature limit and approaches \(1/(2 \sqrt{F^2 \mu^2})\) for \(E \to 0\), which is finite. Thus, \(G_1\) does not exhibit a divergent behavior in the vicinity of \(E = 0\) in the zero-temperature limit.7 For illustration purposes, we compare the functions \(G_1\) and \(G_2\) for various temperatures including the zero-temperature limit in 2.
Let us start the discussion of the flow equation for the previously introduced Silver-Blaze-symmetric version of LPA, see 28 . Its energy functions read \[\label{eq:bosonic-energy-function-sbs-sym} \begin{align} E^{\mathrm{B}}_1 &= \sqrt{k^2 + \tfrac{1}{\bar{\Delta}} \partial_{\bar{\Delta}} U_k(\bar{\Delta}, \mu=0) - F^2 \mu^2} \, , \\ E^{\mathrm{B}}_2 &= \sqrt{k^2 + \partial^2_{\bar{\Delta}} U_k(\bar{\Delta}, \mu=0) - F^2 \mu^2} \,. \end{align}\tag{35}\] To obtain \(U_k(\bar{\Delta}, \mu)\) for this truncation, we first compute \(U_k(\bar{\Delta}, \mu=0)\) at vanishing chemical potential which can then be plugged into the right-hand side of the flow equation 29 for nonzero chemical potential. After that, the resulting equation can simply be integrated with respect to the RG scale \(k\). Provided that the expressions under the square roots of the energy functions 35 are greater than zero, this procedure is in principle well-defined. However, for a given solution \(U_k(\bar{\Delta}, \mu=0)\), positions in \(\bar{\Delta}\)-space may exist (depending on \(\mu\)) for which the expressions under the square roots of the energy functions 35 become negative, especially in the case of SSB, leading to an ill-defined \(U_k(\bar{\Delta}, \mu)\). These points would then need to be excluded from the analysis, spoiling the consistency and applicability of this truncation. Moreover, since the scale-dependent effective potential is not computed self-consistently with this procedure, it may in general not be convex.
These issues can be resolved by using our LPA truncations which take the \(\mu\)-dependence of the scale-dependent effective potential into account. We start by examining the diquark energy functions for the truncation map \({\mathbb{T}}^{\mathrm{LPA}_1}\) which suffers from the double-counting issue as discussed above. In this case, the diquark energy functions read \[\label{eq:bosonic-energy-functions-v1} \begin{align} E^{\mathrm{B}}_1 &= \sqrt{k^2 + \tfrac{1}{\bar{\Delta}} \partial_{\bar{\Delta}} U_k - F^2 \mu^2} \, , \\ E^{\mathrm{B}}_2 &= \sqrt{k^2 + \partial^2_{\bar{\Delta}} U_k- F^2 \mu^2} \,. \end{align}\tag{36}\] These are similar to the ones in 35 . However, \(U_k\) now denotes the full \(\mu\)-dependent scale-dependent effective potential \(U_k(\bar{\Delta}, \mu)\). Plugging these energy functions into the diquark propagator, we find that the presence of the aforementioned singular point in the flow equation results in the following constraints for \(U_k\) in the IR limit: \[\begin{align} \addtocounter{equation}{1}\label{eq:double-counting-2} \partial^2_{\bar{\Delta}} U \geq F^2 \mu^2 \, , \quad \frac{1}{\bar{\Delta}}\partial_{\bar{\Delta}} U \geq F^2 \mu^2 \, , \end{align}\tag{37}\] In words, the double-counting issue prevents the dynamical formation of a nontrivial minimum in the effective potential as associated with SSB. Therefore, from here on, we shall discard this truncation.
Finally, let us discuss our third LPA truncation, i.e., \({\mathbb{T}}^{\mathrm{LPA}_2}\), which does not suffer from the double-counting issue. Here, the diquark energy functions read \[\begin{align} \addtocounter{equation}{1}\label{eq:bosonic-energy-functions-v2} E^{\mathrm{B}}_1 = \sqrt{k^2 + \tfrac{1}{\bar{\Delta}} \partial_{\bar{\Delta}} U_k} \, , \quad E^{\mathrm{B}}_2 = \sqrt{k^2 + \partial^2_{\bar{\Delta}} U_k} \,. \end{align}\tag{38}\] Plugging these energy functions into the diquark propagator, we now find that the singularity of the propagator \(G(E)\) entering the flow equation for the scale-dependent effective potential (see our discussion of 31 above) yields the following constraints for \(U_k\) in the IR limit: \[\begin{align} \addtocounter{equation}{1}\label{eq:convexity} \partial^2_{\bar{\Delta}} U \geq 0 \, , \quad \frac{1}{\bar{\Delta}}\partial_{\bar{\Delta}} U \geq 0\,. \end{align}\tag{39}\] This implies that the scale-dependent effective potential becomes convex as \(k\to 0\), and that a nontrivial ground state as associated with SSB can be dynamically generated in the RG flow. With this at hand, we shall discuss the phase structure of the QDM model in 4.
Let us now discuss the initial condition for the scale-dependent effective potential at the initial RG scale \(\Lambda\). For all discussed versions of LPA, it has the form \[\begin{align} \addtocounter{equation}{1}\label{eq:initial-condition-for-U} U_{\Lambda}(\bar{\Delta}, \mu) = \frac{1}{2}\bar\nu^2 \bar{\Delta}^2 \, , \end{align}\tag{40}\] which is \(\mu\)-independent. Here, we would like to emphasize that the initial condition for the scale-dependent effective average action \(\Gamma_{\Lambda} \equiv S\) should not be confused with the truncation \({\mathbb{T}}\circ \Gamma_\Lambda\) at the cutoff scale \(\Lambda\). Since we discuss the lowest order in the derivative expansion, the truncation \({\mathbb{T}}\circ \Gamma_k\) (even at \(k = \Lambda\)) comes with a kinetic term for the diquark fields, see 22 . This is not the case for the action \(\Gamma_{\Lambda} \equiv S\), see 1 . Thus, we have \(\Gamma_{\Lambda} \neq {\mathbb{T}}\circ \Gamma_\Lambda\) in our present study.8 This discrepancy can be lifted by going to a higher order in the derivative expansion by including the scale dependence of the wave-function renormalization of the diquark fields with initial condition \(Z_{k\to \Lambda}\to 0\).
If we had included a Silver-Blaze-symmetric kinetic term for the diquark fields in the action 1 , which would have corresponded to using \(U_{\Lambda}(\bar{\Delta}, \mu) = \frac{1}{2}(\bar\nu^2 - F^2\mu^2) \bar{\Delta}^2\) as initial condition for \(U_k\) instead of 40 , the squared UV mass \((\bar\nu^2 - F^2\mu^2)\) would decrease for increasing \(\mu\), “boosting" the formation of the diquark condensate. Eventually, in the limit \(F^2\mu^2 \to \bar\nu^2\) corresponding to \(U_{\Lambda}\equiv 0\), the diquark condensate would even diverge, which could, however, be cured by including higher order diquark self-interactions in the initial condition at the UV scale \(\Lambda\). Instead, in our study, the quarks tend to drive the system dynamically into the regime associated with SSB. We emphasize that this is different from models where the bosons are the fundamental degrees of freedom. For example, in purely bosonic models, a negative curvature of the scale-dependent effective potential at its origin is necessarily required at the initial RG scale in order to encounter SSB in the IR limit, see also our discussion of Bose-Einstein condensation in the context of a RBG in 5.1.
Finally, we would like to emphasize that the aforementioned double-counting issue can not be resolved by solely adapting the initial condition. However, it can be resolved by changing the definition of the scale-dependent effective potential in the flow equation. This becomes apparent by introducing a new scale-dependent effective potential \(V_k\) defined as \(V_k = U_k + \tfrac{1}{2} F^2 \mu^2 \bar{\Delta}^2\). Inserting it into the flow equation 29 with the energy functions 38 , we obtain the flow equation with the double-counting issue, where \(U_k\) has been replaced by \(V_k\). Therefore, by initializing the “double-counting flow equation" for \(V_k\) with \(V_\Lambda = \frac{1}{2}\bar\nu^2 \bar{\Delta}^2 + \tfrac{1}{2} F^2 \mu^2 \bar{\Delta}^2\), we obtain the same effective potential \(U\) in the IR limit as that obtained from the flow equation without double counting when we employ the relation \(U_{k = 0} = V_{k = 0} - \tfrac{1}{2} F^2 \mu^2 \bar{\Delta}^2\).
Figure 3: Phase diagram of the QDM for LPA as defined by the truncation \({\mathbb{T}}^{\mathrm{LPA}_2}\). In the left panel, the contour lines are associated with the size of the condensate in LPA (solid lines) and MFA (dashed lines), respectively. The right-hand panel shows a zoom into the first-order region depicted in the left panel. This region opens up at low temperatures and moderate chemical potentials. The red line depicts the first-order phase transition. The red dot is the tricritical point located at \((\mu,T)\approx (188.2\,\mathrm{MeV}, 2.6\,\mathrm{MeV})\). In the shaded area, we find indications for the formation of a scale-dependent effective potential of the type associated with a first-order phase transition but the symmetry is eventually restored in the RG flow, see, e.g., 4 (a)..
Figure 4: Derivative of the scale-dependent effective potential for various RG scales (see left panel) at \(T = 1 \,\mathrm{MeV}\) for \(\mu=190\,\mathrm{MeV}\) (left panel, symmetric phase) and \(\mu=210\,\mathrm{MeV}\) (right panel, phase associated with SSB)..
As discussed in detail in the previous section we shall restrict ourselves to the truncation 22 with the temporal kinetic term for the diquark fields as defined in 25 , i.e., we choose \({\mathbb{T}}^{\mathrm{LPA}_2}\) and compare the results to those from MFA. To this end, we have to specify \(\bar{\nu}^2\) in the action \(\Gamma_\Lambda[\Phi] \equiv S[\Phi]\), see 1 . Inspired by QCD phenomenology, see, e.g., Ref. [61] and former studies on diquark condensation, see Refs. [2], [18], [39], [40], we choose \(\bar{\nu}^2\) such that \(\bar{\Delta}_\mathrm{gs} /\sqrt{2} \approx 100\,\mathrm{MeV}\) at a low temperature of \(T=10\,\mathrm{MeV}\) and \(\mu=350\,\mathrm{MeV}\), where the chiral symmetry is expected to be restored in QCD with two massless quark flavors. For both MFA and LPA, we find \(\bar{\nu}_\mathrm{MFA}^2=0.0585\,\mathrm{GeV}^2\) and \(\bar{\nu}_\mathrm{LPA}^2=0.0575\,\mathrm{GeV}^2\), respectively. A complete list of the (numerical) parameters underlying all figures in the present work can be found in [tab:numerical-data] in 8.
Let us first consider the phase diagrams obtained from MFA and LPA in 3 (a) which are found to be quite similar. Indeed, the contour lines associated with the size of the condensate agree very well at rather large \(\mu\) and small \(T\). While this observation is a trivial consequence of our parameter fixing procedure near \(\mu \approx 350\,\mathrm{MeV}\), the good agreement over a wide range of chemical potentials is still remarkable and states that fluctuation effects are suppressed in this region of the phase diagram, at least with respect to the minimum of the effective potential. By increasing the temperature, we then find that the contour lines of the condensate start to deviate and bosonic fluctuations tend to boost the restoration of the symmetry.
The most interesting difference between MFA and LPA occurs at the phase transition at low temperatures, \(T\lesssim4\,\mathrm{MeV}\). While the phase boundary within MFA is of second-order and bends towards lower \(\mu\) for decreasing \(T\),9 the phase boundary within LPA has a turning point and bends towards higher \(\mu\) for decreasing \(T\) below that point. This can be seen more clearly in 3 (b). From 3 (b), we also deduce that the order of the phase transition changes from second to first order slightly below the turning point, which is reflected by the fact that the contour lines merge at the phase boundary. This implies that a tricritical point occurs at \((\mu,T)\approx (188.2\,\mathrm{MeV}, 2.6\,\mathrm{MeV})\), where first- and second-order phase transition lines merge. Since we are numerically restricted to computations with \(T \gtrsim 1\,\mathrm{MeV}\), see 8, we can not make definitive statements about the shape of the first-order phase transition at even lower temperatures. In any case, this phase structure indicates that the scaling behavior of the condensate 18 in MFA at \(T=0\) is absent when fluctuation effects are taken into account.
The first-order phase transition is generated by an interplay of the BCS singularity10 in the quark sector and the singularity of the diquark propagator functions at low temperatures and moderate chemical potentials. For low temperatures, the strength of the singularity of the diquark propagator functions is reduced (in the sense of the principle of strongest singularity, see Ref. [60]) and simultaneously the BCS singularity controlling the dynamics in the quark sector becomes more pronounced. More precisely, the BCS singularity generates steep negative slopes in \(\partial_{\bar{\Delta}}U_k\) close to \(\bar{\Delta}=0\) for RG scales \(k\lesssim \mu\), see also 1. The singularity in the diquark propagator functions is sensitive to the formation of this negative slope and causes \(\partial_{\bar{\Delta}} U_k\) to exceed zero, resulting in the formation of a small bump in \(\partial_{\bar{\Delta}} U_k\) close to \(\bar{\Delta} = 0\) at \(k \lesssim \mu\). The blue shaded region in 3 (b) shows where such a bump (i.e., a signal of a first-order phase transition) still occurs in the RG flow for nonzero \(k\) in the symmetric phase. The boundary of this region intersects with the phase transition line at the tricritical point. In the symmetric phase, close to the phase transition, the aforementioned bump remains present in form of a plateau in \(\partial_{\bar{\Delta}} U_k\) down to the very small RG scales, see 4 (a), whereas it is fully removed in the broken phase, see 4 (b). More precisely, a comparison of 4 (a) and 4 (b) illustrates that the quark contribution becomes stronger with increasing \(\mu\). This continuously lowers the plateau in \(\partial_{\bar{\Delta}} U\) with increasing \(\mu\), as depicted in 5. Therefore, the diquark condensate \(\bar{\Delta}_\mathrm{gs}\) (as associated with the minimum of \(U_k\) in the IR limit) changes discontinuously at some critical value of \(\mu\) where the height of the plateau becomes zero.
Finally, we would like to mention that a stronger singular behavior of the diquark propagator functions may prevent the formation of a bump associated with a first-order transition in the scale-dependent effective potential. This could potentially change the order of the phase transition in the corresponding regime of the phase diagram from first to second order. However, we expect that a first-order phase transition should also be present for other regulator shape functions in the diquark sector since the Litim regulator already exhibits the strongest singular behavior within the standard class of regulator shape functions in LPA, see Ref. [60].
The considerations presented in 3.3 can also be applied to other models including complex scalar fields at nonzero chemical potential such as the RBG and the QMDM.
In the following, we consider a massive RBG, see, e.g., Refs. [43], [57], [62], [64], [65], with a single complex scalar field \(\varphi\) at nonzero chemical potential and temperature. The super field \(\Phi\) is now given by \(\Phi^{\mathrm{RBG}} = (\varphi, \varphi^*, \boldsymbol{\mu})\). The action reads \[\begin{align} \addtocounter{equation}{1}\label{eq:action-RBG} S[\Phi^{\mathrm{RBG}}] = \int_x \Big(&D_{-}\varphi^* D_{+}\varphi + \vec{\nabla} \varphi^* \vec{\nabla} \varphi \\ &+ m^2 \varphi^* \varphi + \lambda (\varphi^* \varphi)^2\Big) \, , \end{align}\tag{41}\] where the covariant derivatives are \(D_{\pm} = (\partial_0 \pm \boldsymbol{\mu})\), \(m\) is the mass parameter, and \(\lambda\) is the four-boson coupling. The scale-dependent effective potential \(U_k\) is obtained from the projection rule 11 with the homogeneous field configuration \[\begin{align} \addtocounter{equation}{1}\label{eq:scale-dependent-effective-potential-RBS} \Phi^{\mathrm{RBG}}_0 = (&\varphi=\tfrac{1}{\sqrt{2}}\bar{\varphi}, \varphi^*=\tfrac{1}{\sqrt{2}}\bar{\varphi}, \boldsymbol{\mu} = \mu) \, , \end{align}\tag{42}\] where \(\bar{\varphi}\) and \(\mu\) are real numbers. The flow equation for the scale-dependent effective potential reads \[\begin{align} \addtocounter{equation}{1}\label{eq:LPA-flow-equation-RBG} \partial_t U_k &= - v_{d-1} k^{d}\sum_{\chi=\pm1}\Big(1 + \chi\frac{4 \mu^2}{\xi^2_+ - \xi^2_-}\Big) G_2(\xi_\chi)\,. \end{align}\tag{43}\] Here, \(G_2(E)\) is given by 32 with \(F = 1\). The energy functions \(E^{\mathrm{B}}_{1/2}\) entering in \(\xi_{\pm}\) depend again on the variant of LPA used in the concrete calculations, see 3.3. For the same reasons as discussed above for the QDM study, we are solely considering \({\mathbb{T}}^{\mathrm{LPA}_2}\) whose energy functions are given in 38 with \(F = 1\). The flow equation 43 is initialized at \(k = \Lambda\) with the potential \[\begin{align} U_{\Lambda}(\bar{\varphi}, \mu) = \frac{1}{2}(m^2 - \mu^2) \bar{\varphi}^2 + \frac{1}{4}\lambda\bar{\varphi}^4 \, . \end{align}\] Note that this initial potential is \(\mu\)-dependent in contrast to the one chosen for the QDM, see 40 . This \(\mu\)-dependence emerges from the temporal kinetic term for the bosons included in the action \(S\), see 41 and causes \(U_{\Lambda}(\bar{\varphi}, \mu)\) to develop a nonzero minimum already at the scale \(\Lambda\) for sufficiently large \(\mu\). It is precisely this \(\mu\)-dependence of the initial condition that can lead to a nontrivial ground state in the IR limit, even for \(m^2>0\). From a phenomenological standpoint, the presence of this nontrivial ground state is associated with SSB and the formation of a Bose-Einstein condensate in this model, see, e.g., Refs. [57], [66]. Recall that bosonic fluctuations tend to restore the symmetry in the IR limit.
In 6, we show our result for the phase boundary for the parameter choice \(m/\Lambda=0.1\) and \(\lambda=1\). We find \(\mu_{\rm c}/\Lambda=0.1561\), which is the critical value of the chemical potential for Bose-Einstein condensation in the zero-temperature limit. For comparison, we included the results from a previous fRG study [62] a lattice calculation [63], and the mean-field phase boundary [66] given by \(T_\mathrm{MFA}(\mu)=\sqrt{3 (\mu^2-m^2)/\lambda}\) for \(\mu \geq m\). Note that the temperature as well as the chemical potential are given in units of \(\mu_{\rm c}\). In all these studies, the phase transition is found to be of second-order and the corresponding critical temperature increases with the chemical potential. Moreover, we find very good agreement with the previous fRG study. However, our results deviate significantly from MFA results. Notably, our results are in good agreement with those from the lattice calculations at low temperatures, indicating that the truncation map \({\mathbb{T}}^{\mathrm{LPA}_2}\) indeed provides us with a good truncation although it explicitly breaks the Silver-Blaze symmetry. For increasing chemical potential, we then find that our results for the critical temperature start to deviate from the lattice results which may at least partially be traced back to cutoff artifacts.
Finally, we present two examples of RG flows in the Bose-Einstein condensation phase in 7: one corresponding to \(T/\mu_\mathrm{c}=0.064\) and the other one to the zero-temperature limit. Both RG flows break down at (different) nonzero RG scales due to the limitation of the numerical precision in our calculations. For \(T=0\), the RG flow breaks down at higher RG scales than for \(T/\mu_\mathrm{c}=0.064\). This can be traced back to the strength of the singularity of the propagator function \(G_2\) which becomes stronger when the temperature is increased, see our discussion below 34 and also 2. This also implies that the scale-dependent effective potential becomes convex faster in the RG flow which is indicated by the typical flattening of \(\partial_{\bar{\varphi}} U_k\) as illustrated in 7. At the respective (numerical) breakdown scales, the position of the minimum of \(U_k\) (zero of \(\partial_{\bar{\varphi}}U_k\)) is already almost converged in both cases shown in 7. This suggests that our results do not suffer from the numerical instabilities in this respect.
In this subsection we add meson fields \(\phi_i\) (with \(i = 1,2,3,4\)) to our QDM action in 1 which leads us to a QMDM. Models of this type have previously been studied with quarks coming in two colors, see, e.g., Refs. [44], [67]–[69], and three colors, see, e.g., Refs. [24], [27], [69]–[72].
In the following we consider the action \[\begin{align} \addtocounter{equation}{1}\label{eq:action-QMD} S[\Phi^{\mathrm{QMD}}] = \int_x \biggl[&\bar{\psi}_{c} \bigl( {\rm i} \cancel{\partial} - {\rm i} \boldsymbol{\mu} \gamma^0 \bigr) \psi_{c} + \bar{\nu}^2 \Delta_{c}^* \Delta_{c} \nonumber \\ &- h_\Delta \bar{\psi}_{c_1} \gamma^\mathrm{ch} C \tau^2 \Delta_{c_0}^* \frac{{\rm i}}{2} \epsilon_{c_0 c_1 c_2} \bar{\psi}_{c_2} \nonumber \\ &+ h_\Delta \psi_{c_1} C \gamma^\mathrm{ch} \tau^2 \Delta_{c_0} \frac{{\rm i}}{2} \epsilon_{c_0 c_1 c_2} \psi_{c_2} \nonumber\\ &+h_\sigma \bar\psi_c {\rm i} ([\phi_1 + {\rm i} \sum_{i=2}^{4}\tau^{i-1} \phi_i \gamma^5]\big) \psi_c \nonumber\\ &+\tfrac{1}{2} \partial_\nu\phi_i\partial_\nu\phi_i + \frac{1}{2} m^2 \phi_i \phi_i + \frac{1}{4} \lambda (\phi_i \phi_i )^2 \biggr] \; , \end{align}\tag{44}\] where \(\tau^{i}\) with \(i = 1,2,3\) are the three Pauli matrices, \(m\) and \(\lambda\) are the mass and quartic coupling associated with the meson fields, respectively. The Yukawa couplings \(h_\Delta\) and \(h_\sigma\) are associated with the diquark and meson fields, respectively. The super field is defined as \[\begin{align} \addtocounter{equation}{1}\label{eq:super-field-QMD} \Phi^{\mathrm{QMD}} = (\Delta_c, \Delta^*_c, \bar\psi_c, \psi_c, \phi_i, \boldsymbol{\mu}) \, . \end{align}\tag{45}\] A numerical study of this model is beyond the scope of the present work. Here, we shall only show its LPA flow equation as obtained from a straightforward but proper generalization of our truncation map \({\mathbb{T}}^{\mathrm{LPA}_2}\) underlying the QDM study. To this end, we define the projection rule for the scale-dependent effective potential as follows: \[\begin{align} \addtocounter{equation}{1}\label{eq:scale-dependent-effective-potential-QMD} U_k(\sigma, \bar{\Delta}, \mu) = \frac{1}{V_d} \Gamma_k[\Phi=\Phi^{\mathrm{QMD}}_0] \, , \end{align}\tag{46}\] where \[\begin{align} \addtocounter{equation}{1}\label{eq:super-field-evaluation-QMD} \Phi^{\mathrm{QMD}}_0 = (&\Delta_c=\tfrac{1}{\sqrt{2}}\bar{\Delta}\,\delta_{c 3}, \Delta^*_c=\tfrac{1}{\sqrt{2}}\bar{\Delta}\,\delta_{c 3}, \\ &\quad \bar\psi_c=0, \psi_c=0, \phi_i=\sigma\,\delta_{i 1}, \boldsymbol{\mu} = \mu) \, . \end{align}\tag{47}\] Note that the scale-dependent effective potential depends on one more field direction than the QDM model, requiring the use of a suitable framework in numerical studies, see, e.g., Ref. [36].
Applying \({\mathbb{T}}^{\mathrm{LPA}_2}\) to the QMDM, we eventually obtain the flow equation for the scale-dependent effective potential:
\[\begin{align} \addtocounter{equation}{1}\label{eq:LPA-flow-equation-QMDM} \partial_t U_k = &-v_{d-1} k^d \bigg(2(N_\mathrm{c} - 1) \, G_1(E_{\mathrm{G},\Delta}) + (N_\mathrm{f}^2 - 1) \, G_2(E_{\mathrm{G},\sigma}) + G_3(E_{\mathrm{G},\Delta}, E_{\mathrm{R},\Delta}, E_{\mathrm{R},\sigma}, \partial_{\bar{\Delta}} \partial_\sigma U_k)\bigg) \\ &+ v_{d-1} k^d N_\mathrm{f} d_\gamma \bigg(\sum_{\chi = \pm 1} \frac{k}{E^{\mathrm{F}}} \frac{E^{\mathrm{F}}+\chi \mu}{2E^\mathrm{F}_\chi}\bigl(1-2 n_\mathrm{F}(\beta E^\mathrm{F}_\chi)\bigr) + \frac{k}{2E^\mathrm{F}} \Bigl(1- n_\mathrm{F}\big(\beta(E^\mathrm{F}+\mu)\big)- n_\mathrm{F}\big(\beta(E^\mathrm{F}-\mu)\big)\Bigr)\bigg)\, , \end{align}\tag{48}\]
where \[\label{eq:energy-functions-QMDM} \begin{align} E_{\mathrm{G},\Delta}&=\sqrt{k^2 + \frac{\partial_{\bar{\Delta}} U_k}{\bar{\Delta}}}\, , \\ E_{\mathrm{R},\Delta}&=\sqrt{k^2 + \partial^2_{\bar{\Delta}} U_k}\, ,\\ E_{\mathrm{G},\sigma}&=\sqrt{k^2 + \frac{\partial_\sigma U_k}{\sigma}}\, , \\ E_{\mathrm{R},\sigma}&=\sqrt{k^2 + \partial^2_\sigma U_k}\, , \\ E^\mathrm{F}_\chi&=\sqrt{\frac{1}{2} h_\Delta^2 \bar{\Delta}^2 + (E^{\mathrm{F}}+\chi\mu)^2}\, , \\ E^{\mathrm{F}}&=\sqrt{k^2+h_\sigma^2\sigma^2} \, \end{align}\tag{49}\] are the bosonic and quark energy functions and \(G_1(E)\) is given in 31 , \(G_2(E)\) in 32 and \[\begin{align} \addtocounter{equation}{1}\label{eq:new-propa-QMDM} &G_3(E_{\mathrm{G},\Delta}, E_{\mathrm{R},\Delta}, E_{\mathrm{R},\sigma}, \partial_{\bar{\Delta}} \partial_\sigma U_k) \\ & \qquad =k \sum_{j=1}^{3} \frac{\mathcal{N}(r_j)}{\sqrt{r_j} \prod\limits_{l=1, l\neq j}^{3}(r_j-r_l)} \Bigl[ \frac{1}{2} + n_\mathrm{B}(\beta \sqrt{r_j}) \Bigr] \, . \end{align}\tag{50}\] Here, we have introduced the auxiliary function \(\mathcal{N}(r)\), \[\begin{align} \addtocounter{equation}{1} &\mathcal{N}(r) = (-r+E_{\mathrm{R},\sigma}^2)\bigl( (-r+E_{\mathrm{R},\Delta}^2) + (-r+E_{\mathrm{G},\Delta}^2) \bigr) \\ &- (\partial_{\bar{\Delta}} \partial_\sigma U_k)^2 + (-r+E_{\mathrm{R},\Delta}^2) (-r+E_{\mathrm{G},\Delta}^2) + 4 F^2\mu^2 (-r) \, . \end{align}\] In 50 , \(r_j\) are the roots of \[\begin{align} \addtocounter{equation}{1} &\mathcal{D}(r)= - (\partial_{\bar{\Delta}} \partial_\sigma U_k)^2 (-r+E_{\mathrm{G},\Delta}^2) \\ &+(-r+E_{\mathrm{R},\sigma}^2)\bigl( (-r+E_{\mathrm{R},\Delta}^2) (-r+E_{\mathrm{G},\Delta}^2) + 4 F^2\mu^2 (-r) \bigr) \, . \end{align}\] As in the case of the QDM, SSB at large chemical potential is controlled by an interplay of the BCS singularity in the quark sector and bosonic fluctuations. In the IR limit, the singularity in the bosonic sector then ensures convexity of the scale-dependent effective potential.
In the present work, we constructed various truncation maps to define for theories which involve complex scalar fields coupled to a chemical potential. We find that the Silver-Blaze-symmetric version of LPA is not suitable for investigating SSB. One of the Silver-Blaze-violating maps, i.e., the one defined by 22 with 25 , turns out to be best suited for an investigation of SSB in these type of theories. Our discussion of Silver-Blaze-symmetric truncation maps remains a well-defined starting point for a systematic construction of truncations for studies of theories at small chemical potentials at higher order in the derivative expansions.
In 4 and 5.1, we have applied the truncation map 22 to a QDM and to a RBG, respectively. For the QDM, we find a first-order phase transition at low temperatures which becomes a second-order phase transition at a tricritical point. The emergence of a first-order phase transition is the result of a nontrivial competition of the BCS singularity in the quark sector and the singular behavior of the diquark propagator functions at small temperatures. Interestingly, our results indicate the existence of a first-order phase transition from a trivial ground state to a color superconductor at a nonzero critical chemical potential at zero temperature. This suggests that the typical BCS-type exponential scaling behavior 18 of the gap usually observed in MFA studies of this type of models is absent when fluctuations of the order-parameter field are taken into account as done in LPA. For the RBG, we find a second-order phase transition for all considered values of the chemical potential. Notably, our results are found to agree well with those from previous fRG and lattice studies.
Finally, we briefly discussed how a can be constructed from our QDM by suitably including mesonic degrees of freedom and provided a corresponding truncation map, see 5.2. The resulting flow equation for the scale-dependent effective potential may serve as a starting point of future numerical studies of this model which may provide us with a fresh insight into the competition of chiral symmetry breaking and the formation of a color superconductor in dense two-flavor QCD matter.
We thank A. Geißel, P. Guttandin, A. Koenigstein, J. M. Pawlowski, F. Rennecke, B. Schallmo, S. Töpfel, M. Thies, and A. Zorbach for discussions and comments on the manuscript. As members of the fQCD collaboration [73], the authors also would like to thank the members of this collaboration for discussions. This work is supported in part by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) through the Collaborative Research Center CRC 1245 “Nuclei: From Fundamental Interactions to Structure and Stars" – project number 279384907 – SFB 1245, the CRC-TR 211”Strong-interaction matter under extreme conditions" – project number 315477589 – TRR 211, and the State of Hesse within the Research Cluster ELEMENTS (Project No. 500/10.006).
This appendix considers the dependence of the scale-dependent effective potential and thermodynamic observables on the chemical potential in theories with Silver-Blaze symmetry at zero temperature on more general grounds. To this end, we use the notation for the generalized scale-dependent effective potential as introduced in 3.3.
In theories where bosons do not couple to the chemical potential (e.g., the multi-component meson field \(\phi\) in a two-flavor quark-meson model), Silver-Blaze symmetry implies, at the level of the scale-dependent effective average action, that \[\begin{align} \addtocounter{equation}{1} \Gamma_k[\bar{\psi}_c,\psi_c,\phi](\boldsymbol{\mu}) = \Gamma_k[\bar{\psi}_c {\rm e}^{{\rm i} \alpha \cdot}, {\rm e}^{-{\rm i} \alpha \cdot} \psi_c, \phi](\boldsymbol{\mu} - {\rm i}\alpha)\,, \end{align}\] where the symbol \(\cdot\) represents a placeholder for the argument of the transformed fields. Note that we have used 5 and that uncharged scalar fields transform trivially under Silver-Blaze transformations. By evaluating both sides of this relation on a constant meson background field configuration \(\phi=\sigma\) and setting \(\bar{\psi}_c=\psi_c=0\), we arrive at the following relation for the generalized scale-dependent effective potential: \[\begin{align} \addtocounter{equation}{1}\label{eq:Ukrelation} {\boldsymbol{U}}_k(\sigma,\boldsymbol{\mu}) = {\boldsymbol{U}}_k(\sigma,\boldsymbol{\mu} -{\rm i}\alpha)\,, \end{align}\tag{51}\] where we have used the definition \[\begin{align} \addtocounter{equation}{1} {\boldsymbol{U}}_k(\sigma,\boldsymbol{\mu}) = \frac{1}{V_d}\Gamma_k[0,0,\sigma](\boldsymbol{\mu})\,. \end{align}\] When considering the zero-temperature limit, \(\alpha\) is a continuous variable. From 51 we then deduce that \({\boldsymbol{U}}_k\) does not depend on the imaginary part of the chemical potential. Moreover, if the generalized scale-dependent effective potential \({\boldsymbol{U}}_k(\sigma, \boldsymbol{\mu})\) is complex analytic in \(\boldsymbol{\mu}\) for a given tuple \((k, \sigma)\), it must also be independent of the real part of \(\boldsymbol{\mu}\) at that tuple. Looking now at the Wetterich equation, we note that the RG flow of the generalized scale-dependent effective potential is determined by the full scale-dependent two-point function \(\Gamma^{(2)}_k\) evaluated on a homogeneous background field \(\sigma\). In practice, the \(k\)- and \(\sigma\)-dependent poles of this quantity therefore restrict the analyticity of \(\boldsymbol{U}_k\) at most to a strip in the complex plane defined by \(0 \leq \Re(\boldsymbol{\mu}) < \tilde{\mu}(k,\sigma)\) and we will assume the existence of such a strip in the following. Assuming that \(\tilde{\mu}(k,\sigma)\) is a strictly increasing and continuous function in \(\sigma\) for any \(k\), the equation \(\tilde{\mu}=\tilde{\mu}(k,\sigma)\) can be uniquely solved for \(\sigma\) which yields a critical value \(\sigma_{\rm c}(k,\tilde{\mu})\). For \(\sigma > \sigma_{\rm c}(k,\tilde{\mu})\), it then follows that the generalized scale-dependent effective potential is complex analytic in the domain \(0 \leq \Re(\boldsymbol{\mu}) < \tilde{\mu}\) and, together with 51 , it follows that \[\begin{align} \addtocounter{equation}{1}\label{eq:Uk95indep95mu} {\boldsymbol{U}}_{k}(\sigma,\boldsymbol{\mu})= {\boldsymbol{U}}_{k}(\sigma,\boldsymbol{\mu}=0) \, . \end{align}\tag{52}\] Let us now turn to the IR limit \(k \to 0\). Provided that \(\sigma_{\rm gs}({\boldsymbol{\mu}}=0) > \sigma_{\rm c}(k=0,\tilde{\mu})\), where \(\sigma_{\rm gs}({\boldsymbol{\mu}}=0)\) is the position of the ground state of the generalized effective potential in the vacuum, 52 implies that \(\sigma_{\rm gs}({\boldsymbol{\mu}})\) does not depend on \(\boldsymbol{\mu}\) in the corresponding chemical potential domain. Conversely, the equation \(\mu_c=\tilde{\mu}(k\!=\! 0,\sigma_{\rm gs}({\boldsymbol{\mu}}\!=\!0))\) defines a critical value of the chemical potential below which the ground state and thus also the pressure \(p=-{\boldsymbol{U}}_{k\to 0}(\sigma_\mathrm{gs}(\boldsymbol{\mu}),{\boldsymbol{\mu}})\) as well as all thermodynamic observables derived from it remain independent of \(\boldsymbol{\mu}\) (which is dubbed Silver-Blaze phenomenon/property [46] in the case of real \(\boldsymbol{\mu}\)).11
For example, this is a well-known zero-temperature feature of models in which bosons do not couple to the chemical potential, such as the two-flavor quark-meson model and the Gross-Neveu model in the zero-temperature limit.12
For a theory where bosons couple to the chemical potential (as it is case for the QDM and the RBG), the situation is different. On the level of the scale-dependent effective average action, Silver-Blaze symmetry implies in this case that \[\begin{align} \addtocounter{equation}{1} & \Gamma_k[\bar{\psi}_c,\psi_c,\Delta_c^{\ast},\Delta_c](\boldsymbol{\mu}) \\ & \quad = \Gamma_k[ \bar{\psi}_c {\rm e}^{{\rm i} \alpha \cdot}, {\rm e}^{-{\rm i} \alpha \cdot} \psi_c, \Delta^*_c {\rm e}^{-{\rm i} F \alpha \cdot}, {\rm e}^{{\rm i} F \alpha \cdot} \Delta_c](\boldsymbol{\mu} - {\rm i}\alpha)\,, \end{align}\] where we have again used 5 . However, by evaluating now both sides on a constant diquark background field configuration \(\Delta^{\ast}_c=\Delta_c=\bar{\Delta}/\sqrt{2}\) and setting \(\bar{\psi}_c=\psi_c=0\), we arrive at \[\begin{align} \addtocounter{equation}{1}\label{eq:Ukcharged} & {\boldsymbol{U}}_k(\bar{\Delta},\boldsymbol{\mu}) \\ & \quad = \frac{1}{V_d} \Gamma_k[0, 0, (\bar{\Delta}/\sqrt{2}) {\rm e}^{-{\rm i} F \alpha \cdot}, {\rm e}^{{\rm i} F \alpha \cdot} (\bar{\Delta}/\sqrt{2}) ](\boldsymbol{\mu} - {\rm i}\alpha) \end{align}\tag{53}\] with \[\begin{align} \addtocounter{equation}{1} {\boldsymbol{U}}_k(\bar{\Delta},\boldsymbol{\mu}) = \frac{1}{V_d} \Gamma_k[0, 0, \bar{\Delta}/\sqrt{2} , \bar{\Delta}/\sqrt{2} ](\boldsymbol{\mu})\,. \end{align}\] Relation 53 corresponds to relation 51 , which is derived for a theory where the bosons do not couple to the chemical potential. In contradistinction to theories of the latter type, however, it is not possible to draw any conclusions from 53 regarding the dependence of the generalized scale-dependent effective potential on the chemical potential for our theory with charged bosons. In this case, the generalized scale-dependent effective potential may therefore carry a dependence on the chemical potential at any RG scale \(k\), even if the underlying scale-dependent effective average action is invariant under Silver-Blaze transformations. For example, this is the case for the generalized scale-dependent effective potential of the QDM in MFA, see Eq. 17 . In particular, the generalized scale-dependent effective potential may exhibit a nontrivial dependence on \(\text{Im}({\boldsymbol{\mu}})\). However, in the absence of SSB,13 where the ground state is associated with \(\bar{\Delta}_\mathrm{gs}=0\), we can deduce from 53 that \[\begin{align} \addtocounter{equation}{1} {\boldsymbol{U}}_k(0,\boldsymbol{\mu}) = {\boldsymbol{U}}_k(0,\boldsymbol{\mu} -{\rm i}\alpha)\,. \end{align}\] Thus, in the absence of SSB at \(T=0\), \({\boldsymbol{U}}_k(0,\boldsymbol{\mu})\) does not depend on the imaginary part of the chemical potential. If, additionally, \({\boldsymbol{U}}_k(0,\boldsymbol{\mu})\) is complex analytic in \(\boldsymbol{\mu}\) for some \(k\), then \({\boldsymbol{U}}_k(0,\boldsymbol{\mu})\) is even independent of the real part of the chemical potential. If this holds in the limit \(k\to 0\) for \(0\leq \Re({\boldsymbol{\mu}}) < \mu_{\rm c}\), where \(\mu_{\rm c}\) again denotes the extent of the complex analytic domain in the direction of \(\Re({\boldsymbol{\mu}})\), we conclude that thermodynamic observables in the absence of SSB are independent of the chemical potential below a critical value \(\mu_c\).
Let us finally comment on the existence of a nonzero critical value \(\mu_c\) of the chemical potential. In our discussion above, we assumed complex analyticity in the domain \(\Re(\boldsymbol{\mu})< \mu_{\rm c}\). From a thermodynamic standpoint, it is reasonable to assume the existence of such a domain. Indeed, the chemical potential is the change in free energy when a charged particle (i.e., a particle that couples to a given chemical potential, such as a quark or diquark in case of a quark chemical potential) is added to or removed from the system. Therefore, the density of the charged particles can become finite only for \(\mu \geq \mu_{\rm c}\), where the critical value \(\mu_{\rm c}\) is set by the vacuum pole mass of the lightest charged particle in the spectrum of the theory,14 see Ref. [45] for a detailed discussion. For \(\mu < \mu_{\rm c}\), the density remains zero and is independent of the chemical potential. Since the density is directly related to the first derivative of the grand canonical partition function with respect to chemical potential, the latter (and therefore also the pressure) does not depend on the chemical potential for \(\mu < \mu_{\rm c}\). As discussed above, this invariance of the grand canonical partition function under a shift of the chemical potential below \(\mu_{\rm c}\) can be understood as a consequence of the Silver-Blaze symmetry of the effective action and complex analyticity of the ground-state value of the generalized effective potential in \(\boldsymbol{\mu}\).
The flow equation 29 is a partial differential equation that can be recast into a continuity equation by performing a field-derivative on both sides, see Refs. [33], [37]. In this continuity formulation, the MFA contribution is a source term, the pure \(E^{\mathrm{B}}_1\)-dependent bosonic contribution is an advection term and the mixed \(E^{\mathrm{B}}_1\) and \(E^{\mathrm{B}}_2\)-dependent bosonic contribution is an advection-diffusion term. This formulation enables the use of the well-established KT scheme [75] to solve the flow equation 29 . Therefore, we always worked on the level of the field derivative of the scale-dependent effective potential rather than the potential itself. As a side effect, this formulation establishes a well-defined boundary condition \(\partial_{\bar{\Delta}} U_k(\bar{\Delta}=0, \mu)=0\) at \(\bar{\Delta}=0\) for \(T>0\). At low \(T\), this is where the “relic” of the BCS singularity determines the slope of \(\partial_{\bar{\Delta}} U_k\) and where it is therefore crucial to control the dynamics very well.
\begin{ruledtabular}
\setlength\extrarowheight{8pt}
\begin{tabular}{c c c c c c c}
Figure & Resolution & Range & \gls{UV} parameter & $k_{\mathrm{IR}}$ & $tol$ &\\
\hline
\ref{fig:MF-BCS} & $0.1$ & $300$ & $\bar{\nu}^2=58.5\times10^3$ & $5$ & $10^{-14}$ &\\
\ref{fig:LPA-V1-phase-diagram-contour-full} & $2$ & $2\times10^3$ & $\bar{\nu}_\mathrm{LPA}^2=57.5\times10^3$ & $75$ & $10^{-10}$ &\\
\ref{fig:LPA-V1-phase-diagram-contour-full} & $2$ & $2\times10^3$ & $\bar{\nu}_\mathrm{MFA}^2=58.5\times10^3$ & $75$ & $10^{-10}$ &\\
\ref{fig:LPA-V1-phase-diagram-contour-first-order-regime} & $1$ & $10^3$ & $\bar{\nu}^2=57.5\times10^3$ & $13$ & $10^{-10}$ &\\
\ref{fig:LPA-V1-first-order-flow} & $1$ & $10^3$ & $\bar{\nu}^2=57.5\times10^3$ & $14$ & $10^{-14}$ &\\
\ref{fig:LPA-second-order} & $1$ & $10^3$ & $\bar{\nu}^2=57.5\times10^3$ & $14$ & $10^{-14}$ &\\
\ref{fig:BEC} & $0.5$ & $10^3$ & $m^2=10^4, \lambda=1$ & $1$ & $10^{-10}$ & \\
\ref{fig:RBG-flow} & $0.5$ & $10^3$ & $m^2=10^4, \lambda=1$ & $14$ & $10^{-12}$
\end{tabular}
\end{ruledtabular}
\caption{\label{tab:numerical-data}List of model and numerical parameters for all plots presented in this work.
All dimensionful variables are given in powers of $\MeV$ according to their mass dimension.
The \gls{UV} scale has always been chosen to be $\Lambda=1000\MeV$ and $rtol=atol=tol$.
``Resolution" and ``Range" refer to the grid in field space.
To resolve the phase boundary (``zero height contour line") in \gls{LPA} more accurately, see \ref{fig:LPA-V1-phase-diagram-contour-full}, we employed a range of $10^3$ and a resolution of $1$ in field space with $k_\mathrm{IR}=20$.}
We implemented the KT scheme in a semi-discrete way, i.e., we discretized the field space, but used a continuous RG scale, which leads to a coupled system of ordinary differential equations [74]. We treated the mixed
advection-diffusion term as advective in the sense of the KT scheme. We then used Python’s solve_ivp [76] in Python 3 [77] (using various libraries [76], [78], [79]) to solve
this system with “LSODA” as numerical time stepper with the parameters \(rtol\) and \(atol\) for its relative and absolute error. Some of the results have been cross-checked with
Mathematica
[80]. A plain list of model parameters and numerical parameters for all figures shown in this work can be found in
[tab:numerical-data]. Note that in some domains associated with SSB, as shown in 3 (a), but far from the phase boundary, we stopped the RG flow at comparatively large values of \(k_{\text{IR}}\). This is justified because the RG flow close to the minimum of the scale-dependent effective potential “freezes out” early in the RG flow anyhow, due
to strong SSB in these domains.
However, dealing with mixed advection-diffusion terms is challenging also for the KT scheme. In regions of the phase diagram, where the quark contribution creates a deep well, i.e., at low temperatures and large chemical potentials, we find that numerical noise occurs when we consider very high resolutions. More precisely, the solution \(\partial_{\bar{\Delta}} U_k(\bar{\Delta}, \mu)\) develops spurious oscillations around the minimum of this well, i.e., where the mixed contribution changes from being diffusion dominated to advection dominated. However, the influence of these artifacts on the results presented in our present work appear to be negligible. In fact, the results appear to be properly converged coming from lower grid resolutions.
Note that there are applications in which phenomenologically meaningful fields, such as the \(\omega^0\) meson, enter a theory in the same way as the imaginary part of a complex-valued chemical potential see, e.g., Ref. [42].↩︎
Here, we also assume that the generalized truncation reduces to the original truncation at real \(\boldsymbol{\mu} = \mu\).↩︎
Note that the potential \(U_k\) itself has an \(\text{O}(2N_\mathrm{c})\) symmetry at \(k=\Lambda\), see 1 , which is preserved during the RG flow, because the right-hand side of 14 solely depends on the \(\text{O}(2N_\mathrm{c})\) invariant \(\Delta_{c}^* \Delta_{c}\).↩︎
Since the nonanalytic behavior is encoded in Eq. 17 in terms which become independent of the regulator in the limit \(k^\prime \to 0\), it follows that an expansion of the scale-dependent effective potential about \(\bar{\Delta}=0\) is ill-defined for \(k^\prime \leq \mu\) for any regulator from the class of standard \(3d\) regulators. This issue may be circumvented by integrating out fluctuations about the Fermi surface but this then requires the use of regulators that break the Silver-Blaze symmetry by construction [45]. Note that, as discussed in Ref. [55], an interchange of the integration of the flow equation over the RG scale \(k\) and an expansion in powers of a background field is in general delicate. However, this does not underlie the appearance of a divergent curvature of the effective potential at the origin in the present case.↩︎
Recall the transformation of the chemical potential in Eq. 5 and the fact that the generalized scale-dependent effective potential is in general built up from powers of the chemical potential up to arbitrarily high orders.↩︎
Note that alternative generalizations of the scale-dependent effective potential with full \(\mu\)-dependence (e.g., involving only \(\Re(\boldsymbol{\mu}\))) may also be invariant under Silver-Blaze transformations. However, unlike 19 , they explicitly violate complex analyticity and consequently cannot be utilized to study features such as the Silver-Blaze property.↩︎
In contrast to neglecting the contributions associated with the function \(G_1\) in the RG flow, neglecting the terms associated with \(G_2\) would imply that the RG flow loses the aforementioned “self-healing” property [60] in the limit \(T\to 0\), i.e., the system may potentially flow into a point with \(E=0\) for \(k>0\).↩︎
However, we are still fulfilling the property 9 for the only projection occurring in our truncation, i.e., the scale-dependent effective potential.↩︎
In MFA, we expect that the phase boundary tends to zero asymptotically for \(\mu\to 0\), since SSB occurs at all \(\mu\) for \(T=0\) for our choice of model parameters, see 3.2.↩︎
In the zero-temperature limit, the BCS singularity manifests itself in a nonanalytic behavior of the scale-dependent effective potential associated with a singular behavior of its curvature at the origin. This singular behavior is removed at nonzero temperature. However, at sufficiently low temperatures, the physics is still dominated by the “relic" of this singularity present in the zero-temperature limit. Here and in the following, we shall for simplicity also refer to this aspect as BCS singularity.↩︎
The Silver-Blaze phenomenon/property should not be confused with the Silver-Blaze symmetry of the effective action. The Silver-Blaze phenomenon/property refers to the fact that observables do not depend on the chemical potential below some critical value. This phenomenon necessarily requires the Silver-Blaze symmetry of the effective action. However, the symmetry of the effective action under Silver-Blaze transformations does not imply the existence of a finite Silver-Blaze region where thermodynamic observables are independent of the chemical potential.↩︎
For instance, in the Gross-Neveu model using a spatial Litim regulator in MFA, \(\boldsymbol{U}_k(\sigma, \boldsymbol{\mu})\) remains \(\boldsymbol{\mu}\)-independent, provided that \(k^2 + h_\sigma^2 \sigma^2 > \Re(\boldsymbol{\mu})^2\), see, e.g., Refs. [55], [74]. The resulting critical value is \(\mu_{\mathrm{c}} = h_\sigma \sigma_{\mathrm{gs}}(\boldsymbol{\mu}=0)\). However, note that \(\boldsymbol{U}_{k \to 0}(\sigma, \boldsymbol{\mu})\) may still exhibit a \(\boldsymbol{\mu}\)-dependence for \(\sigma < \sigma_{\mathrm{gs}}(\boldsymbol{\mu}=0)\), even within \(0 \leq \Re(\boldsymbol{\mu}) < \mu_c\).↩︎
We also tacitly assume the absence of any explicit symmetry breaking in terms of, e.g., external sources.↩︎
Strictly speaking, in the relativistic systems considered in this work, it is not the density but the difference in the densities of charged particles and their antiparticles.↩︎