Analyticity and the phase diagram of QCD

Some consequences of the analyticity of the free energy (pressure) of QCD at finite chemical potential are deduced. These include a method for numerical exploration of the full phase diagram by a novel use of simulations at imaginary chemical potential to extract Yang-Lee zeroes of the grand-canonical partition function. We make use of, and comment on, CPT symmetries, positivity of non-linear susceptibilities and the finiteness of screening lengths. We also comment on the structure of zeroes expected for the usual picture of the phases of QCD, following a discussion of the physics of imaginary chemical potential.

Motivated by the complex phase structure of QCD at finite chemical potential, µ, predicted by effective theories and perturbative analyses [1], there has been a spate of work extending lattice computations to finite µ [2,3,4,5,6]. Most of these methods extrapolate data obtained at µ = 0 out to finite µ using, implicitly or explicitly, a Taylor expansion in µ of the free energy. As a result, these expansions cannot be continued beyond the nearest phase transition to µ = 0. Much of the interesting phase structure of QCD then lies beyond their scope. This paper explores consequences of the analyticity of the free energy (guaranteed by time reversal invariance) as a function of complex µ to suggest a numerical technique for exploring the QCD phase diagram. It turns out that analyticity also has many other consequences, some of which are touched upon.
The plan of this paper is the following. First it is shown how analyticity is guaranteed by time reversal symmetry. This material is known to practitioners, but serves to introduce notation and several concepts I use later. Next I write down the canonical product representation for the partition function and make contact with the Yang-Lee theory of phase transitions. A numerical technique is developed to explore the phase diagram using simulations at imaginary chemical potential, which goes beyond using Taylor series for analytic continuation. The need for a different approach is discussed next, with remarks on the phase diagram.
The partition function of QCD, i.e., SU(N c ) gauge theory, with N f flavours of fermions, each subjected to a common chemical potential, µ, is where M is the Dirac operator, S g the gauge part of the action and the temperature T enters the action through boundary conditions on the fields. The free energy, F and the pressure are defined by At any purely real or imaginary value of the chemical potential, the free energy (or the pressure) is purely real [7]. Is it possible then for the free energy to be an analytic function of a general complex µ? Recall that this requires the real and imaginary parts of F , F r and F i respectively, to satisfy the Cauchy-Riemann equations, implying that they are conjugate harmonic functions. This rules out F i = 0 unless F is constant.
The answer is that it is entirely possible to have a harmonic function (i.e., a function satisfying Laplace's equation) which is zero along the real and imaginary axes and is nontrivial elsewhere. Since all even powers of the complex number µ satisfy these requirements, it seems that analyticity of F (µ) (at fixed T ) requires that it be even in µ.
Now recall that the transformation µ → −µ can be compensated by time-reversal, i.e., by interchanging particles and anti-particles. However, this definition is arbitrary; in the absence of CP violating terms in the action, thermodynamics remains unchanged by such a relabelling. The physics is therefore unchanged by changing the sign of µ. In other words, CP symmetry (same as time reversal, by the CPT theorem) implies that F (µ) is an even function of µ, i.e., only even powers of µ appear in the Taylor expansion. Since the whole chain of logic above can be inverted, this means that CP symmetry is necessary and sufficient for the analyticity of F , and its reality along the coordinate axes of µ in the complex plane.
The simplest use of analyticity is to expand the pressure in a Taylor series [5,6] and to determine the coefficients by direct lattice simulations at µ = 0 [6,8] or by fitting to results obtained at imaginary µ [4]. However, this approach is limited, since it only gives information on the phase connected to µ = 0. The first phase transition encountered stops the extrapolation, and more machinery is required to get around this barrier. This is what we develop in this paper.
At real µ, the path integral measure in eq. (1) is complex [9] and importance sampling in a Monte Carlo procedure fails although the integral is real. From the argument above, it seems that by grouping together configurations related by CP, it should be possible to construct a real weight by summing over each such "CP orbit". Since CP symmetry forms a Z 2 group, each orbit consists of exactly two configurations, and the sum of the weights in these two configurations (for real µ) must be real. The simplest example is the Gibbs model The determinant is- where w = qU and ma/2 = sinh m ′ a/2. The lattice spacing is a = 1/T N t . The partition function is which is real, although the weight is complex. In this simple model a CP transformation induces the mapping U → U † . Then, summing over each CP orbit one gets which is real. In the partition function one should compensate the double counting by dividing by a factor of two. This is immaterial for expectation values.
In a more realistic model, the gauge configurations connected by a CP transformation are harder to construct. It turns out to be easier to implement the symmetry transformation C (which is the same as PT) which is local and maps each link matrix U → U † . The T part of the transformation can be used as before to prove that the free energy is even in µ.
Moreover, this symmetry transformation preserves the values of S g and Det M (see eq. 1) separately. As a result, for the Taylor expansion of Det M summed over each PT orbit- one has ℜ a i = ℑ a i = 0 for all odd i. In other words, summed over PT orbits, the measure is even in µ.
Since this summation over a symmetry orbit leads to a weight that is the sum of two determinants, it is clear that the measure cannot be rewritten using local pseudo-fermion fields. Disappointingly, these symmetries turn out to be unexploitable in this form for molecular dynamics algorithms paralleling those which are used at zero chemical potential [11]. However, by ensuring analyticity, they lead to other numerical approaches to the problem, as I show next.
From the Hamiltonian expression for the partition function, (whereĤ is the Hamiltonian andN the number operator) it is clear that Z has a periodicity of 2πT in the imaginary part of the chemical potential- Functions with such periodicity are best analysed in terms of the fugacity, z = exp(µ/T ).
Entire functions in µ are regular in the complex z plane without the origin, i.e., the punctured z plane [12]. Then the partition function can be written in the canonical product form [13] Z where z n (T ) are the zeroes of Z, and F is the regular part of the free energy. CP symmetry constrains the set of zeroes, {z n } to be symmetric under inversions in the unit circle. Branch points in the full free energy, F = −T log Z, then do not appear in the regular part F . The periodicity of Z in µ also implies that The appropriate setting for further analysis is the Yang-Lee theory of phase transitions [14].
According to this theory, phase transitions occur at points in the parameter space where dense sets of zeroes develop into pinch points. In the present understanding of the phase diagram of QCD, one expects pinch points on the real µ axis at small T where interesting new phases of QCD appear [1]. We shall come back to this after a digression on ideal gases which throws more light on the question of periodicity in the imaginary part of µ.
An ideal gas has the free energy- which does not have the required periodicity; no polynomial in µ can. By examining the Hamiltonian formulation of an ideal gas, we can show where the periodicity is hidden. The free energy can be written as Since z is unchanged under µ → µ + 2πiT , this has explicitly the required periodicity.
By the substitution y = exp x, the integral can be written in terms of fourth order Euler Polylogarithm functions [15]- which also explicitly retains the periodicity. In fact, each of the terms has this periodicity explicitly. This property is lost in a Taylor expansion of the polylogs in µ/T - where the remaining terms are odd in the expansion parameter. The reason is that in taking the logarithm of the argument one makes a choice of the Riemann sheet, which is precisely equivalent to working within a single strip of width 2πiT . When the sum in eq. (14) is expanded in a Taylor series in µ, it is clear that every term cancels, apart from those which appear in eq. (12). Periodicity can be restored by averaging this expression over all Riemann sheets of the logarithm, the result of which is expressed concisely in eq. (14) [16]. In general, interactions modify the integral expression in eq. (13), without spoiling the periodicity. This upsets the conspiracy which cancels all higher terms in F , leading to non-vanishing values for general Taylor coefficients, i.e., the non-linear number susceptibilities [6].
We return now to considerations of phase transitions in QCD through the Yang-Lee mechanism. The set of zeroes, z n (T ), of the partition function contains the full information needed to identify all the phase transitions in the theory through the Yang-Lee mechanism [14]. But for all this to make sense, one needs to check that a thermodynamic limit exists, and that, in this limit, the pressure is a continuous and monotonically increasing function of z.
This is Theorem I of [14], and is proven there in the non-relativistic limit when interparticle potentials obey two conditions-first that there should be a limit to the number of particles that can be accommodated in a finite box, and second that the interactions between them be of finite range. Non-intuitive phenomena can occur when these conditions are violated [17]. All work on QCD assumes the validity this theorem. However, it is useful to first examine the basis for this assumption.
In the low temperature phase, where the carriers of baryon number are the baryons themselves, the problem is essentially non-relativistic, since the baryon mass is much greater than the temperature. There is ample experimental evidence that the inter-nucleon potential satisfies both conditions necessary for the theorem to hold-the baryon-baryon interaction has a hard core repulsion and a short ranged Yukawa interaction at long distances. In the hightemperature phase the problem is relativistic, since one expects that the carriers of baryon number are the quarks, whose masses are much less than the temperature. Quarks and antiquarks are freely created and destroyed, and, as is well-known, the correct reformulation of the problem is to work in ensembles where the excess (or deficit) of quarks over antiquarks is fixed (the canonical ensemble) or is conjugate to the chemical potential (the grand canonical ensemble). The second requirement, of finite range of interactions, is straightforward.
The gauge interaction, which could be long-ranged, is known to be screened in the electric sector [18]. In the magnetic sector there is no evidence for screening, but there is evidence of confinement and hence a mass-gap [19]. All other screening lengths are smaller in the continuum limit [8]. While these phenomena are best established for N c = 2 and 3, they are expected to hold also for other N c .
If we accept these arguments, or otherwise directly assume the first theorem of [14]-that the pressure is a continuous increasing function for z ≥ 1, then an immediate consequence is that the Taylor expansion coefficients of P (T, z) in z at fixed T are non-negative. Now, the generalised susceptibilities [6], are Taylor coefficients of P (T, µ) in an expansion in µ.
Since the expansion in µ can be obtained from that in z = exp(µ/T ) by expanding out the exponential, which is itself a convex function, it follows that the sum of all generalised susceptibilities of a given order are strictly positive.
This theorem of Yang and Lee is expected to hold for a continuum theory in the thermodynamic limit of infinite volume. It is also expected to hold for a continuum theory in a finite volume much bigger than the range of interactions, even when the limit is taken with a sequence of cutoffs, holding the physical volume fixed. On the other hand, at any intermediate step in this procedure, the convexity argument may fail if some of the screening masses are far from their continuum limits and not sufficiently small. Negative values of some of the susceptibilities may then be observed at finite volume. However, convexity is regained either by increasing the spatial volume or on taking the continuum limit [20].
where Z(T, iφ) denotes the grand canonical partition function in an imaginary chemical potential µ = iφ. One can invert this formula to write where N s is the number of components of the Dirac spinor (N s = 4 in four dimensions).
Note that on any lattice with V sites on each spatial slice, Z is a polynomial in z of degree where the coefficients are real and we can write |α N − 1| = ǫ N , |α 0 − 1| = ǫ 0 and |α i | = ǫ i (for 0 < i < N). We model the errors ǫ i as being independent random numbers drawn from Gaussians of width σ. For sufficiently large statistics σ is small, and the deviations of the roots of eq. (18) from z n , δ n , can be treated as linear in the ǫ i . It is then a straightforward exercise to show that δ n = 0 and the RMS error, δ 2 n − δ n 2 = σz n /N. Thus the problem is statistically well conditioned. The estimates of zeroes are unbiased and the errors fall inversely with the square root of the statistics and inversely with the system volume.
We turn to the question of phase transitions at purely imaginary µ, since this yields further insight into the locations of the zeroes of the partition function. We note at the outset that crucial convexity theorems (such as the second law) fail at imaginary µ, and hence thermodynamics in its usual sense does not apply. The idea of [3] is to use the peak in a response function, χ, (i.e., second derivative of the free energy with respect to an intensive parameter) to identify putative critical points at imaginary µ and finite volume, analytically continue them to real µ, and then take the infinite volume limit. Note several subtleties in the argument. First, for a function χ(µ 2 ), the extrema for imaginary chemical potential, µ 2 * < 0, are obtained by solving for the zeroes of the first derivative. This does not necessarily give the extrema at real chemical potential, since the solution of χ ′ (µ 2 * ) = 0 does not guarantee that χ ′ (−µ 2 * ) also vanishes. A second, related subtlety can be seen in a double Taylor expansion for any response function at finite volume- where the Taylor coefficients c nm are independent of T and µ [3]. At the critical end-point, T E (µ E ) this series must sum to the divergent quantity |1 − T /T E (µ E )| −γ (with γ > 0) in the infinite volume limit. There may also be a critical point, (T ′ E , µ ′ E ) at an imaginary chemical potential, (µ ′ E ) 2 < 0. Clearly the sum in eq. (19) differs at the two points, i.e., even if χ(T E /T c , µ E /T E ) and [24]. In the context of our earlier arguments, it is more natural to look for accumulation points of zeroes of the partition function for imaginary chemical potential. This leads to the third subtlety, and a genuine difference between real and imaginary chemical potential due to the lack of convexity mentioned earlier. The theorem of Yang and Lee which we have earlier quoted allowed us to express the partition function as a finite (at all finite volumes) Laurent series in z with non-negative coefficients (eq. 17). This constraint on the coefficients ensures convexity and prevents zeroes at real positive z, i.e., for real µ. Only in the limit of infinite volume can the zeroes pinch the real line and cause a phase transition. However, this same structure allows, (and, in some models, may force) the occurrence of zeroes on the unit circle in z (i.e., for imaginary µ) even at finite volume.
Examples abound. We point out the extension of the Gibbs model to SU(N c ) gauge group where the partition function has been computed [25], and found to be This is of the form in eq. (17) The zeroes on the imaginary µ axis give rise to logarithmic branch point singularities of the free energy evaluated at imaginary chemical potential. The quark number density, all quark number susceptibilities, the energy density and all its derivatives diverge at all these points for all values of N c . These divergences are thermodynamically meaningless since they occur without the necessity of taking a thermodynamic or continuum limit. At finite quark mass, the zeroes remain on the unit circle but move about as long as Q(T, 0) ≤ 1. When the mass becomes larger, these zeroes leave the imaginary axis, but they are still visible as maxima of these observables (see Figure 1), and can be interpreted as crossover points. The effect of mass seems to be one of the main differences between the perturbative and strong coupling computations of [7].
These crossovers are generic at imaginary chemical potential-Potts models (or O(N) models) with an ordering field in one direction and an imaginary transverse field show such behaviour, as do massive fermions subjected to imaginary chemical potential. Removing the first ordering field converts them into branch point singularities. Both pose dangers to numerical simulations, since these non-thermodynamical singularities are easily misinterpreted as strong first order transitions, or developing critical points. If one indeed misidentifies such points, then finite size scaling studies are also misleading, since one tries to apply inappropriate homogeneity relations arising from thermodynamic considerations. The idea of locating the zeroes of the partition function using data collected at imaginary chemical potential is to bypass these problems altogether.
The connection of the results of [7] and the phase diagram of QCD at real µ seems to be rather subtle. Here we make a hypothesis which is consistent with the observed absence of a critical point with Z 3 symmetry [26] and with present understanding of the physics of the Z Nc phases [27]. There are also two insets showing the positions of Yang-Lee zeroes which might lead to such phase diagrams. In this figure only the pinch points are constrained-the remainder of the lines 'a' and 'b' are allowed to twist, merge or intersect without changing the phase diagram [29]. The important point is that in a range of temperatures, there are two pinch points on the real axis [30]. As the temperature rises, the pinch point on the right, labelled 'b', moves to higher and higher µ r , whereas the one on the left, labelled 'a', lifts away from the real axis at the critical end point, T E . Its effect can still be felt at the crossover temperature, T c , seen in finite temperature simulations at µ = 0. Above T c the line 'b' continues to move outwards, but there remain isolated N c zeroes on or near the unit circle (drawn with a dotted line), depending on the quark mass, which give rise to the crossovers (F and G) observable at imaginary µ.
We end by summarizing the contents of this paper. The complex analytic structure of the free energy (pressure) of QCD was investigated. It was shown that CP symmetry makes the free energy an analytic function of the chemical potential µ at all temperatures T . This symmetry is unexploitable to give an algorithm for a Monte Carlo procedure that works through a molecular dynamics algorithm. Exploiting analyticity, it is possible to write a canonical product representation of the partition function, and made contact with the Yang-Lee theory. It was argued that analyticity, monotonicity and continuity of the pressure implies that non-linear susceptibilities of arbitrary order must exist and be nonnegative at all T . It was further argued that this is connected with electric screening and magnetic confinement in high temperature QCD. I showed how to use canonical partition functions, derived from Monte Carlo simulations at purely imaginary chemical potential, to obtain information on the phase structure of QCD. Finally, a plausible phase diagram of QCD was used to anticipate the behaviour of the zeroes that such numerical techniques might give rise to.