Universal aspects of the phase diagram of QCD with heavy quarks

The flavor dependence of the QCD phase diagram presents universal properties in the heavy quark limit. For the wide class of models where the quarks are treated at the one-loop level, we show, for arbitrary chemical potential, that the flavor dependence of the critical quark masses-for which the confinement-deconfinement transition is second order-is insensitive to the details of the (confining) gluon dynamics and that the critical temperature is constant along the corresponding critical line. We illustrate this with explicit results in various such one-loop models studied in the literature: effective matrix models for the Polyakov loop, the Curci-Ferrari model, and a recently proposed Gribov-Zwanziger-type model. We further observe that the predictions which follow from this one-loop universality property are well satisfied by different calculations beyond one-loop order, including lattice simulations. For degenerate quarks, we propose a simple universal law for the flavor dependence of the critical mass, satisfied by all approaches.

Understanding the properties of strongly interacting matter in extreme conditions of temperature, density, magnetic field, etc., is a question of topical interest with various implications, e.g., in early Universe cosmology and astrophysics [1]. This relates to some of the most profound aspects of quantum chromodynamics (QCD), namely, the physics of (de)confinement and of chiral symmetry breaking/restoration. Past, present, and upcoming heavy-ion collision experiments at RHIC, CERN, FAIR, NICA and J-PARC, accompanied by an intense theoretical effort, aim at unravelling the phase diagram and the thermodynamic properties of QCD [2,3].
First-principle calculations based on numerical lattice simulations have demonstrated a crossover between a mostly confined to a mostly deconfined phase as the temperature T is increased [4]. Such calculations are, however, restricted to small baryon chemical potential µ B due to the infamous sign problem that prevents the use of standard Monte Carlo algorithms [5]. A major open question of both theoretical and experimental research programs concerns the possible existence of a first order transition line in the (T, µ B ) plane, ending at a critical point [6,7]. Continuum approaches can avoid the strong sign problem of the lattice but necessarily rely on approximations or model building that must be tested in situations where lattice results are under control.
In this context, it is of interest to study QCD-like theories, by varying parameters such as the number of colors and flavors, the quark masses, etc. Investigating the phase structure of the theory in such an extended parameter space, beyond its own theoretical interest, may contribute to our understanding of the phase diagram of the physical theory. An example is the investigation of the Columbia plot, which aims at describing the phase structure of the SU(3) color theory with two degenerate quark flavours u and d and one strange quark s, as a function of the quark masses M u = M d and M s [8][9][10].
The most commonly accepted scenario at zero chemical potential and in the limit of vanishing quark masses is that of a first order transition governed by the restoration of chiral symmetry at large temperatures [8]. This transition weakens as one increases the quark masses and it eventually becomes second order for critical values of the latter. The physical point lies beyond this critical line in the Columbia plot, where the transition is a crossover. A question is, therefore, to study how this critical line evolves with nonzero chemical potential [9]. Also of interest is the limit of infinite quark masses, corresponding to the pure Yang-Mills theory. For three colors, the latter presents a first order confinementdeconfinement phase transition in temperature [11]. Again, this transition weakens for large but finite quark masses and eventually becomes second order for critical values of the quark masses. Contrarily to the light quark case, the heavy quark system can be simulated also at nonzero chemical potential because the sign problem can be evaded by using a large mass expansion [12]. In that case, the critical line is known to shrink toward larger quark masses as the chemical potential increases.
Not only does QCD with heavy quarks provide a framework where one can explicitly study the role of the chemical potential, but it also yields a stringent test case for approximate continuum approaches based either on first principle calculations or on models of the QCD dynamics. Existing calculations include truncations of the Dyson-Schwinger equations (DSE) or functional renormalisation group techniques [13], one-and two-loop perturbative calculations in the Curci-Ferrari (CF) model (a massive extension of the Faddeev-Popov (FP) Lagrangian in the Landau gauge) [14,15], and effective models for the Polyakov loop, the order parameter of the confinementdeconfinement transition [7,16,17].
In the present Letter, we discuss generic features of the phase diagram of QCD with heavy quarks. We show arXiv:1805.10015v1 [hep-th] 25 May 2018 that, for a wide class of continuum approaches, where the quarks are treated at the one-loop level, various aspects of the critical line are, in fact, independent of the gluon dynamics. This is, for instance, the case for the flavor dependence of the critical quark masses. Also the defining equation of the corresponding critical line (surface) in the Columbia plot is independent of the chemical potential. Finally, the critical temperature and the critical values of the order parameters are constant along the critical line, as observed in lattice simulations for the former [12].
We illustrate these features on various results from the literature, such as the CF model at one-loop [14] and the matrix model of Ref. [17]. We also provide results for a model of the glue dynamics recently proposed in Refs. [18,19], based on the Gribov-Zwanziger (GZ) quantization scheme, where one restricts the path integral over gauge fields to the first Gribov region, for which the FP operator is positive definite [20]. We also show that our generic findings for the critical line extend beyond oneloop as exemplified by the results of the CF model at two-loop order [15] and by the DSE results of Ref. [13]. Finally, we show that the universality mentioned above implies a simple and robust scaling law for the flavor dependence of the quark masses on the critical line that is well satisfied by all existing continuum approaches and also by the lattice results. This allows us to predict the critical quark masses for increasing number of flavorsprovided confinement is not lost -in all these approaches.
Let us first consider the case of a vanishing quark chemical potential µ = 0 and N f degenerate heavy quark flavors of mass M for simplicity. Descriptions of the QCD dynamics where quarks are included at the one-loop level yield the potential is the modified Bessel function of the second kind and β −1 is the temperature. The Polyakov loop is related to the free energy F q of a static quark in the thermal bath of gluons as ∝ exp(−βF q ) and vanishes identically in the confined phase, where F q = ∞. Here, we assume that the potential v glue is confining, that is, it admits a minimum at = 0 at zero temperature. For each value of N f , the critical values for , β and M on the upper boundary line of the Columbia plot are determined by solving the system of equations The critical values for and β are determined from the two equations in (3) We deduce, in particular, that, even though the critical values c , β c and R N f depend on the considered model for v glue , the relation between the values of R N f for various values of N f is model-independent, at one-loop order.
The above discussion easily generalises to the case of nondegenerate quark masses, where one concludes that the equation determining the critical surface in the space of quark masses is universal, i.e., independent of v glue . For instance, the critical line in the plane (M u = M d , M s ) of the Columbia plot for 2 + 1 (heavy) flavors is determined by the equation where only the number on the right hand side, which determines the location of the critical line, depends on v glue , through Eq. (2). We can use, alternatively, 3f . This trivially generalises to more flavors.
The same strategy applies to the case of nonzero quark chemical potential µ = µ B /3. The effective potential now depends separately on the Polyakov loops and¯ corresponding to quarks and anti-quarks, respectively, and reads with a = ∂ 2 V | c and b = ∂ ∂¯ V | c . The first two equations rewrite As in the previous case, the last two equations in (8), involving two or more and¯ derivatives, only concern v glue and define two functions (β) and¯ (β) independent of both N f and µ. Together with these, the ratio of the two equations (9), determines an N f -independent value for β c (µ) -or, more directly, µ c (β) -and in turn N f -independent critical values c (µ) = (β c (µ)) and¯ c (µ) =¯ (β c (µ)). As before, the only source of flavor dependence comes from the prefactor of the function f in Eq. (9). Taking the ratio of any of these for two different values of N f , we conclude that Eq. (4) still holds, for any value of µ. Again, this can be generalised to the case of nondegenerate quark masses and we conclude, following the same lines of argumentation, that the critical values c (µ), c (µ), and β c (µ) are constant along the whole critical surface and that the flavor dependence of the latter is insensitive to v glue and to the value of µ. For instance, in the 2+1 flavor case, the critical line is defined by Eq. (5). Only the actual location of the critical line depends on v glue and on µ, through the determination of R 3 on the right-hand side of Eq. (5). This is demonstrated in Fig. 1.
Finally, Eq. (9) gives the µ-dependence of R N f in terms of the functions β c (µ), c (µ) and¯ c (µ), defined above. We also emphasize that the above discussion is valid for both real and imaginary chemical potential. In this respect, it is interesting to mention that the phase diagram at imaginary chemical potential features a tricritical point, which results in a particular scaling of both β c and R N f when approaching βµ = iπ/3 [22]. The latter can be analytically continued to real µ and is actually well reproduced by lattice data and by various continuum approaches. The existence of a tricritical point is intimately related to the Roberge-Weiss symmetry of the potential for βµ = iπ/3, which reads V ( ,¯ ) = V (e i2π/3¯ , e −i2π/3 ). To see this, it is convenient to change to the real variables x and y defined by = e iπ/3 (x − iy) and¯ = e −iπ/3 (x + iy), in terms of which the Roberge-Weiss symmetry rewrites y → −y. Defining x(y) from the condition ∂ x V = 0 and evaluating the potential along this trajectory, one obtains a reduced potential V r (y) with Z 2 symmetry y → −y. The tricritical point corresponds to the cancellation of both the second and the fourth derivatives of V r (y) at y = 0. Below, we determine the tricritical point both by following the boundary line of the Columbia plot as βµ approaches iπ/3 or directly from the reduced potential at βµ = iπ/3.
It is interesting to test the above considerations on various one-loop approaches available in the literature. These include the matrix model of Ref. [17], the CF model at one-loop order [14], as well as the GZ type model of Refs. [18,19]. In this latter case, no results are available for the Columbia plot. We produce them here after briefly reviewing the approach.
The model for the glue potential proposed in Refs. [18,19] is inspired from the GZ quantization in the vacuum, where the restriction of the path integral to the first Gribov region is implemented by means of a suitably adjusted Gribov parameter. The extension of this approach to finite temperature proposed in [18,19]  The function W has been computed to one-loop order in these references, where it was expressed in terms of a background gauge field βĀ µ (x) = δ µ0 r j t j , with t j the generators of the Cartan subalgebra of the gauge group. For SU(3), these are the two diagonal Gell-Mann matrices λ (3) /2 and λ (8) /2. At this order of approximation, the background field is related to the Polyakov loops as The color modes κ correspond either to one of the Cartan directions j or to any of the roots α of the gauge algebra. In Ref. [19], three scenarios for the color dependence of the Gribov parameters have been proposed, all compatible with the symmetries of the problem: GZ1: γ κ = γ , ∀κ . GZ2: γ κ = γ 1 if κ = j and γ α = γ 2 , ∀α . GZ3: γ κ = γ 1 if κ = j and γ α = γ −α .
We stress that, when computing the derivatives of the glue potential in Eqs. (7) and (8), one needs to take into account the -and¯ -dependence of the Gribov parameters. Taking derivatives of the gap equation ∂ γ W = 0, one gets, for example in the case GZ1, ∂ γ = −∂ ∂ γ W/∂ 2 γ W . Other derivatives are obtained in a similar way. This easily generalises to the scenarios GZ2 and GZ3.
We compute the values of R N f for various N f in each scenario and compare them to lattice results. At the order of approximation considered here, these quantities do not depend on the value of the Gribov parameter at zero temperature and are, therefore, a stringent test for the various scenarios. Our results are gathered in Table I and compared to other approaches. We observe that, among the various one-loop approaches considered here, the degenerate GZ approach (GZ1) gives the best results as compared to the lattice values. We mention that the background dependence of the Gribov parameters is crucial in obtaining such values. Neglecting the latter yields the values refered to as GZ0 in Table I. We also note that the non-degenerate case (GZ3) is completely discarded : For N f = 1, we find a first order phase transition, irrespectively of the value of the quark mass and the values for R 2 and R 3 are way too small. One easily checks that the universal equation (4) is accurately satisfied, given the precision, for all one-loop results in Table I. In contrast, the lattice values do not satisfy the scaling (4) well. As we discuss below, they follow a somewhat different scaling, where the function f (x) is replaced by a simple exponential.
In Table I, we also quote the results from two calculations beyond one-loop order, namely, the two-loop perturbative calculation in the CF model of Ref. [15] and the DSE results of Ref. [13]. In this case, the direct comparison of critical masses is complicated by nontrivial mass renormalisation beyond leading order. As argued in Ref. [15], this can be partially absorbed in the ratios R N f /R N f . In particular, these ratios for DSE agree well with the other continuum approaches despite the low values of R N f . Note that the same argumentation does not save the scenario GZ3, firstly, because there is no reason to expect such mass renormalisation effect at one-loop order and, secondly, because the ratio R 3 /R 2 is bad anyway as compared to all other continuum approaches.
In the same vein, as discussed in Ref. [15], the lattice also brings additive mass renormalisations because it explicitly breaks chiral symmetry and more suitable quantities to compare are ratios of differences, e.g., In Table I, we quote the results for Y 3 , which is surprisingly stable for all the approaches (continuum and lattice). As far as continuum models of the form (1) are concerned, this can be understood as a simple consequence of the scaling law (4). First, since R N f is always large (except in the GZ3 case but this scenario is discarded), the Bessel function can be approximated by its asymptotic behavior and one easily checks that up to relative corrections of O(R −2 1 ). The above result yields Y 3 ≈ ln 3/ ln 2 ≈ 1.58, in excellent agreement with our results in Table I.
It is remarkable that the results for Y 3 for all approaches beyond one-loop agree well with the prediction (13). We can understand the lattice value from the N fscaling observed in Ref. [23] valid for all values of µ. This yields the scaling law exp(R N f − R N f ) = N f /N f , from which one immediately obtains Eq. (13) using the fact that the lattice value of β c is found insensitive to N f . We conclude that the scaling law (13) is quite robust and universal (independent of the actual gluon dynamics) and can thus be used to predict the values of R N f for N f ≥ 3, given R 2 and R 1 in the various approaches. We have checked that these predictions describe well the actual results in the one-loop models GZ1, GZ2, and CF1, as expected. We have also verified that the same holds for the two-loop values of R N f in CF2. The results for an imaginary chemical potential, and more precisely along the tricritical line at βµ = iπ/3, are gathered in Table II. We verify that the law (4) and its consequences as described above are well verified by the values of R N f for one-loop models. Again, the GZ1 scenario seems to be giving the best values to date, although they are in this case pretty close to the ones from the matrix model of Ref. [17].
As anticipated, the values of Y 3 from all approaches agree well with the prediction (13). This is to be expected from the above discussion for both one-loop models and Lattice results. The DSE value is strickingly good despite the low values of R N f and the not so good values of the ratios R N f /R N f . This supports the conjecture of the universal character of Y N f . Acknowledgement: We thank Daniel Kroff for related