Locate QCD Critical End Point in a Continuum Model Study

With a modified chemical potential dependent effective model for the gluon propagator, we try to locate the critical end point (CEP) of strongly interacting matter in the framework of Dyson-Schwinger equations (DSE). Beyond the chiral limit, we find that Nambu solution and Wigner solution could coexist in some area. Using the CornwallJackiw-Tomboulis (CJT) effective action, we show that these two phases are connected by a first order phase transition. We then locate CEP as the end point of the first order phase transition line. Meanwhile, based on CJT effective action, we give a direct calculation for the chiral susceptibility and thereby study the crossover.


Introduction
Quantum chromodynamics (QCD) has two important features, namely, dynamical chiral symmetry breaking (DχSB) and color confinement. It is generally believed that with increasing temperature or baryon number density strongly interacting matter will undergo a phase transition from the hadronic matter to the quark-gluon plasma (QGP) which is expected to appear in the ultrarelativistic heavy ion collisions. These two phases are generally referred to as the Nambu-Goldstone phase and the Wigner phase concerning its realization of chiral symmetry. A central goal of the worldwide program in relativistic heavy ion collisions is to chart the phase diagram of QCD in the plane of nonzero temperature (T ) and chemical potential (µ) and to locate the position of the critical end point (CEP) [1][2][3].
In the Nambu-Goldstone phase the quarks are confined and chiral symmetry is dynamically broken which leads to a large dynamical quark mass. In the Wigner phase which is thought to be connected with QGP, chiral symmetry is partially restored and quarks are not confined. Theoretically, these two phases are described by two different solutions, the Nambu-Goldstone solution and the Wigner solution of the dressed quark propagator. The existence and properties of these two solutions in the chiral limit (i.e., the current quark mass m = 0) have been widely investigated in the framework of Dyson-Schwinger equations (DSE) approach of QCD. For example, in Ref. [4] the authors studied this problem by taking the infrared part of the gluon propagator in the Maris-Tandy model [5] and inputting it directly into DSE at finite temperature and chemical potential. It is found that in the chiral limit there is an area on the T − µ plane on which both the Nambu and Wigner solutions coexist. However, there remain two questions. Firstly, in the real world, the current quark masses (m) for the u and d quarks are small but nonzero. Actually, many studies have shown that the absence of current quark mass would exhibits a second order phase transition instead of a crossover in certain areas of QCD phase diagram, which results in giving a tri-critical end point(TEP) instead of a CEP [6,7]. Here, to be closer to the real world, we would take m = 5 MeV. Secondly, in earlier work employing the Maris-Tandy model one simply ignores the µ-dependence of the gluon propagator. However, when µ is large, the influence of µ to the gluon propagator should be important. Recently, authors of Refs. [8,9] have noticed this problem and modified the gluon propagator. They have added terms in their own way to suppress the influence of the gluon propagator at large µ. It is reasonable since it characterizes the weakening of interaction between quarks as µ goes up. Using such a modified gluon propagator, the authors of Ref. [10] have found multi-solution of the quark gap equation at T = 0 and µ = 0 beyond the chiral limit (Similarly, the authors of Refs. [11][12][13] studied the multi-solutions of the quark gap equation with nonzero current quark mass in the case of finite temperature and chemical potential by introducing some modification to the normal Nambu-Jona-Lasinio model). The main purpose of this paper is to generalize the study in Ref. [10] to the case of nonzero T and µ and thereby locate the position of CEP by means of DSE.
In this work, we only consider two degenerate light flavors u and d, and hence µ = µ u = µ d ≈ 1 3 µ B throughout this paper.

Gap equation and its solutions
In Euclidean space and under rainbow-ladder approximation,the quark gap equation at finite T and µ reads whereω n =(2n + 1)πT + iµ.
Here we have put the regularization mass scale at infinity so that all renormalization constants are 1. This is allowed since the model used by us renders all the integrations convergent. Here we note that G −1 ( p,ω n ) also depends on T and µ. It can be expressed as [10,14]: Substituting Eq. (2.2) into Eq. (2.1), one obtains the coupled integral equations for complex scalar functions A, B and C. Here we need to specify the form of the model gluon propagator g 2 D µν (Q) with Q = ( p − q,ω n −ω l ). In this work we adopt the model gluon propagator proposed in [8,10] This model is based on the Maris-Tandy model, which is successful in describing the properties of mesons when used with the rainbow-ladder approximation, especially for the ground state pseudo-scalar and vector mesons [5,15]. The parameters D 0 and ω can be fixed by the masses and electro-weak decay constants of the pions and ρ mesons. However, it was found that for ω between [0.3, 0.5] GeV, D 0 which satisfies ωD 0 = (0.8GeV) 3 gave almost the same results [16]. Later in this paper, we will fix ω further. The inclusion of e −αµ 2 /ω 2 does not affect the result at µ = 0. The parameter α characterizes the weakening of interaction between quarks with increasing µ. In the perturbative calculations of finite-temperature QCD, the running coupling constant decreases with µ. So we could expect a decrease in the coupling in the nonperturbative region. Thus α controls the rate of decrease. In this paper, we employ this form of suppression to investigate the CEP.
To locate the CEP, we mainly investigate the phase transition near the CEP. This brings convenience in numerical computation: when T > 120M eV , we only need to calculate a few Matsubara frequencies. Here, we take 40 Matsubara frequencies into consideration and obtain stable results. The parameters are chosen as ω = 0.44 and α = 0.4, which will be explained later. The solutions are plotted in Fig. 1.
From the behavior of the B function, it can be seen that partial restoration of chiral symmetry takes place in two ways. At large T , when µ goes up, the Nambu phase solution goes smoothly into Wigner phase, whereas at small T , two solutions could coexist at the same µ. We will investigate thermodynamic quantities with these two solutions in the next section.

Critical End Point
In this section, we will draw the QCD phase diagram corresponding to Fig. 1. We use the existence of zero pressure difference between the two phases as the criterion for first order phase transition. CEP is located as the end point of the first order phase transition line. Later, we will show the results of chiral susceptibility as a supplement. Quark condensate is also studied.

CJT effective Action and First Order Phase Transition
Cornwall-Jackiw-Tomboulis effective action [17] has been used at finite temperature and density to calculate the partition function [14,18]. Here, we would like to emphasize that it requires that the gluon propagator should not depend on the quark propagator explicitly. Since we are using the rainbow-ladder truncation and a model for the gluon propagator, the CJT effective action is valid here. Actually, in this case the solution of the gap equation is the stationary point for the CJT effective action [19]. So our calculation of the gap equation and the pressure are consistent with each other in this framework. The pressure difference between the two phases is where the subscript N and W stand for N ambu and W igner, respectively. Using the CJT effective action, we have where the trace operation is over color, flavor, Dirac and coordinate indices all through this paper. By Fourier transformation, we can rewrite Eq. (3.2) in the momentum space along with P W (T, µ) . Then Eq. (3.1) becomes 3) The first order phase transition takes place when ∆P(T, µ) = 0. The evolution of ∆P(T, µ) with fixed T and increasing µ is plotted in Fig. 2. Picking out all the points (T, µ) in the T − µ plane where ∆P = 0, we can find the first order phase transition line, which is shown in Fig. 3. As we can see, in this case CEP is located at (T E , µ E ) = (0.142 GeV, 0.095 GeV) .
Here we would like to point out, as can be seen from Fig. 2, while the area for the existence of multi-solution gets narrower with increasing T , there is always one point where ∆P(µ) = 0. The area finally vanishes at the CEP. This fact justifies our use of the criterion for locating CEP due to its consistency with the general phase diagram.

Chiral susceptibility and Crossover
When the chemical potential is lowered to less than µ E , heating the system would make the Nambu solution change continuously into Wigner solution. In this case, we can use various kinds of susceptibility to characterize what is happening here. One kind of susceptibility, i.e., the chiral susceptibility, is often used in lattice QCD [3], Nambu-Jona-Lasinio(NJL) model [20,21], Dyson-Schwinger equation [22,23] and other approaches.
Once the CJT effective action is known, one can in principle derive any thermodynamical quantity from it. In this paper, we will use the CJT effective action to derive directly the chiral condensate and the corresponding chiral susceptibility. The chiral condensate is considered to be the order parameter for the chiral phase transition, which is defined as:   where ψ ψ = ψ ψ u + ψ ψ d . Substituting Eq. (3.2) into Eq. (3.4), and using ∂ ∂α (detX) = detX · Tr(X −1 X ∂α ), we obtain It should be noted that this is a direct definition and works both in the chiral limit and beyond the chiral limit, so we can calculate the quark condensate in our case (m = 5 MeV) directly. Now it is interesting to compare this definition with the generally used one Here, in order to eliminate the ultraviolet divergence, one has subtracted a term G 0 by hand (see for example, Ref. [24]), whereas in Eq. (3.4) the chiral condensate has been obtained directly from the CJT effective action, where this subtraction term appears automatically. Actually, since we have obtained the effective action in Eq. (3.2), it is easier to take the partial derivative directly. Note that the effective action alone is divergent. What is convergent is the difference of two properly chosen effective actions. For our purpose here,  4) is then straightforward. The result is shown in Fig. 5. As we can see, when µ < µ E the quark condensate decreases continuously as T goes up. This can be compared with the case given in Ref. [10] where a first order phase transition take place, i.e., the quark condensate drops discontinuously at µ c .
The chiral susceptibility can be treated similarly. It is defined as Taking the second order partial derivative of the CJT effective action, we can obtain Fig. 5. It can be seen that when µ < µ E , the chiral susceptibility undergoes a smooth change with increasing T . This is where crossover takes place. When approaching the CEP, the curve gets sharper. The peak tends to go to infinity near the CEP, which is a second order phase transition point. We pick out the coordinates of the maximums and obtain the dashed line in Fig. 3.

Conclusions
Taking m = 5 MeV, we are left with two independent parameters: ω from the Maris-Tandy model and α which controls the decrease of coupling with respect to increase of µ. Note that α does not affect the results when µ = 0. Lattice QCD has given many results in this case [26][27][28]. To locate the CEP in our work, we adopt the result that T c = 172±4 MeV. As can be seen from Fig. 6, in our case ω = 0.44 reproduces T c = 168 MeV. After ω is fixed, we try to find the appropriate value for α. Here, we refer to the µ c at zero temperature given by various works [1,29,30]. The values of α and the relevant µ c can be seen from Fig. 7. This result is expected in advance since dynamical chiral symmetry breaking happens at a relatively large coupling. Table. 1 shows our results for different values of α. At first sight it is surprising that µ E goes up as α increases. This is quite different from the case of µ c . However, if one takes a closer look, one will find a decrease in T E , which actually compensates for the increase of µ E . Hence, the total effect is that increasing α would make the CEP rotate clockwise.

Summary
In this paper we generalize the study in Ref. [10] to investigate the QCD phase diagram at finite temperature and density. The gap equation is solved beyond chiral limit and an area for the existence of multi-solution is found. We then use CJT effective action which is consistent with our gap equation to calculate all the thermodynamical quantities. Although the CJT effective action alone is divergent, we could obtain convergent results for pressure difference, quark condensate and chiral susceptibility. The first order phase transition and crossover are both studied, and from this the CEP is identified and located.