Lattice study of thermodynamic properties of dense QC 2 D
N. Astrakhantsev, V. V. Braguta, E.-M. Ilgenfritz, A. Yu. Kotov, A. A. Nikolaev
LLattice study of thermodynamic properties of dense QC D N. Astrakhantsev,
1, 2, ∗ V. V. Braguta,
2, 3, 4, 5, † E.-M. Ilgenfritz, ‡ A. Yu. Kotov,
2, 4, 5, § and A. A. Nikolaev ¶ Physik-Institut, Universit¨at Z¨urich, Winterthurerstrasse 190, CH-8057 Z¨urich, Switzerland Institute for Theoretical and Experimental Physics NRC “Kurchatov Institute”, Moscow, 117218 Russia Moscow Institute of Physics and Technology, Dolgoprudny, 141700 Russia Bogoliubov Laboratory of Theoretical Physics, Joint Institute for Nuclear Research, Dubna, 141980 Russia National University of Science and Technology MISIS, Leninsky Prospect 4, Moscow, 119049 Russia Department of Physics, College of Science, Swansea University, Swansea SA2 8PP, United Kingdom
In this paper we study thermodynamic properties of dense cold SU (2) QCD within lattice simu-lation with dynamical rooted staggered quarks which in the continuum limit correspond to N f = 2quark flavours. We calculate baryon density, renormalized chiral and diquark condensates for var-ious baryon chemical potentials in the region µ ∈ (0 , µ ∈ (0 , µ >
540 MeVthe system becomes sufficiently dense and ChPT is no longer applicable to describe lattice data.For chemical potentials µ >
900 MeV we observe formation of the Fermi sphere, and the system issimilar to the one described by the Bardeen-Cooper-Schrieffer theory where the the diquarks play arole of Cooper pairs. In order to study how nonzero baryon density influences the gluon backgroundwe calculate chromoelectric and chromomagnetic fields, as well as the topological susceptibility. Wefind that the chromoelectric field and the topological susceptibility decrease, whereas the chromo-magnetic field increases with rising of baryon chemical potential. Finally we study the equation ofstate of dense two-color quark matter.
PACS numbers: 12.38.Gc, 12.38.Aw
I. INTRODUCTION
Study of Quantum Chromodynamics at finite baryondensity is an important research topic of modern physicswhich is closely connected to various problems in astro-physics and cosmology. Experimental studies of QCD atfinite baryon density can be carried out in heavy ion col-lision experiments. In particular, the region of the phasediagram with high temperature and small baryon den-sity is well explored at the Large Hadron Collider (LHC)and Relativistic Heavy Ion Collider (RHIC), while thephysical programs of future Facility for Antiproton andIon Research (FAIR) and Nuclotron-based Ion ColliderFacility (NICA) are focused on large baryon density andsmall temperature.At the moment, the theoretical understanding of theQCD phase diagram in ( µ, T ) plane is rather schematic,since the most powerful approach, lattice simulation ofQCD, cannot be directly applied in the region of finitedensity due to the sign problem [1]. Numerous latticeattempts to overcome the sign problem provide reliableinformation only in the region of small baryon density [2].In the absence of straightforward results from latticesimulation of QCD, one applies different analytical ap-proaches to study the ( µ, T ) phase diagram: mean fieldapproaches [3–9], the method of Dyson-Schwinger equa-tions and the renormalization group [10–13], the large– N c ∗ [email protected] † [email protected] ‡ [email protected] § [email protected] ¶ [email protected] approach [14], perturbative QCD [15, 16] and others. Al-though the results obtained theoretically are important,it is rather difficult to estimate the reliability of thesepredictions.One of the possible ways to shed light on the propertiesof dense media is to apply lattice simulation to theorieswhich are similar to QCD but are not plagued by thesign problem. Although such theories differ from realQCD in some aspects, it is believed that these QCD-liketheories can provide important information common fordense media in general. The most popular choices arethe QCD at finite isospin density [17–21] and the two-color QCD at finite baryon density [22, 23]. This paperis devoted to lattice study of the dense two-color QCD.The two-color QCD at finite chemical potential hasbeen studied with lattice simulations quite intensively,see, e.g. [24–35] and references therein. Mostly, these pa-pers are aiming at the study of the phase diagram of two-color QCD in the region of small and moderate baryondensities.The phase structure of dense two-color QCD and itsproperties were studied in our previous papers [36–39],where lattice simulations were carried out at a relativelysmall lattice spacing a = 0 .
044 fm. In this paper we alsoemploy this spacing. Compared to previous studies atlarger lattice spacings, this allows us to extend the rangeof accessible values of the baryon density, up to quarkchemical potential µ > a r X i v : . [ h e p - l a t ] J u l study how nonzero baryon density influences the proper-ties of the gluon background we calculate chromoelectric,chromomagnetic fields and the topological susceptibility.In addition, we shall study the equation of state of densetwo-color QCD.The manuscript is organized as follows. In the next sec-tion we describe our lattice set-up and details of the cal-culation of the observables under consideration. In sec-tion III we present our results on fermionic observables.In section IV we study how nonzero baryon density mod-ifies the properties of the gluon background. Section IVis devoted to our study of the equation of state of densetwo-color QCD. Finally, in the last section we discuss ourresults and draw the conclusion. II. DETAILS OF THE CALCULATION
In this section we briefly describe our lattice setup.More details can be found in papers [29, 37]. In ourlattice study we used the tree level improved Symanzikgauge action [40, 41]. For the fermionic degrees of free-dom we used staggered fermions with an action of theform S F = (cid:88) x,y ¯ ψ x M ( µ, m ) x,y ψ y + λ (cid:88) x (cid:0) ψ Tx τ ψ x + ¯ ψ x τ ¯ ψ Tx (cid:1) , (1)where ¯ ψ , ψ are staggered fermion fields, M ( µ, m ) x,y is thestandard Dirac operator for staggered fermions, m is thebare quark mass. The chemical potential µ is introducedinto the Dirac operator through the multiplication of thelinks along and opposite to the temporal direction byfactors e ± µa , respectively.In addition to the standard staggered fermion actionwe add a diquark source term [24] to equation (1). Thediquark source term explicitly violates U V (1) and allowsto observe diquark condensation even on finite lattices,because this term effectively chooses one vacuum fromthe family of U V (1)-symmetric vacua.In the present study we are going to investigate a the-ory with the partition function Z = (cid:90) DU e − S G · (cid:0) det( M † M + λ ) (cid:1) , (2)where S G is the tree level improved Symanzik gauge ac-tion. In the continuum, the action (2) corresponds to N f = 2 dynamical fermions.The results presented in this paper have been ob-tained in lattice simulations performed on a 32 lat-tice for the set of chemical potential values in the re-gion aµ ∈ (0 , . am = 0 . m π =741(15) MeV ( am π = 0 . ± . r = 0 . µ = 0 amounts to √ σ = 476(5) MeV at a = 0 .
044 fm. Typically, one carries out numerical simulations at afew nonzero values of λ and then extrapolates to zero λ .Notice, however, that numerical simulations in the regionof large baryon density are numerically very expensive.For this reason, in this paper we have chosen a differentstrategy. Most of our lattice simulations are conductedat a single fixed value λ = 0 . λ/m (cid:39) .
1. How-ever, in order to check the λ -dependence of our resultsfor chemical potentials aµ = 0 . , . , . , . , . λ = 0 . , . λ/m (cid:39) . λ → λ parameter (see discussion below).In our paper we are going to calculate the followingfermionic observables • The diquark condensate: a (cid:104) qq (cid:105) = 1 N s N t ∂ (log Z ) ∂λ , (3) • The chiral condensate: a (cid:104) ¯ qq (cid:105) = a (cid:104) ¯ q iα q iα (cid:105) = 1 N s N t ∂ (log Z ) ∂ ( ma ) ; (4) • The quark number density: a n q = 1 N s N t ∂ (log Z ) ∂ ( µa ) ; (5)The baryon density is a conserved quantity and it doesnot require renormalization. The chiral and diquark con-densates require renormalization. To this end, we are go-ing to follow the renormalization procedure analagous to µa . . . . . . . . . n q / T µ, MeV
Figure 1. Quark number density as a function of chemical po-tential. Dashed red line represents the fit by (24), for detaileddiscussion see the paragraph after eq. (26). µa . . . . . . . . Σ ¯ qq µ, MeV
Figure 2. Renormalized chiral condensate (6) as a function ofchemical potential. Dashed red line represents the fit by (23),for detailed discussion see the paragraph after eq. (26). the paper [20] where the authors studied QCD at nonzeroisospin density :Σ ¯ qq = m m π F (cid:2) (cid:104) ¯ qq (cid:105) µ − (cid:104) ¯ qq (cid:105) (cid:3) + 1 (6)Σ qq = m m π F (cid:2) (cid:104) qq (cid:105) µ − (cid:104) qq (cid:105) (cid:3) (7)To get rid of the additive divergences in the chiraland diquark condensates we subtract the correspondingobservable (cid:104) ¯ qq (cid:105) and (cid:104) qq (cid:105) at zero µ = 0. Since thequark mass is renormalized multiplicatively, m r = Zm ,multiplicative divergence falls out in the combination m r ∂/∂m r = m∂/∂m and the quantity m ( (cid:104) ¯ qq (cid:105) − (cid:104) ¯ qq (cid:105) )has a well-defined continuum limit. Other factors in theformula (6) are introduced in such way, that the chiralcondensate Σ ¯ qq is 1 at µ = 0 and 0 if the chiral symmetryis fully restored. To do it we also use the Gell-Mann–Oaks–Renner relation: m π = m (cid:104) ¯ qq (cid:105) / F [23], where F is the constant in front of the kinetic term of the ChiralPerturbation Theory, in the leading order the pion decayconstant f π = F/ µ = 0 thestaggered Dirac operator M ( µ = 0 , m ) has the follow-ing form M ( µ = 0 , m ) x,y = Q x,y + maδ x,y , where theoperator Q satisfies Q † = − Q . The partition func-tion (2) at zero µ = 0 takes the following form: Z = (cid:82) DU e − S G · (cid:0) det( Q † Q + ( ma ) + λ ) (cid:1) and depends onthe quark mass and λ only via the factor m a + λ .Thus, the diquark source parameter λ and the quark If instead of the pion condensate and the pionic source termin QCD at finite isospin density [20] one considers the diquarkcondensate and the diquark source term in dense two-color QCD,both theories look similar in their properties. µa . . . . . Σ qq µ, MeV
Figure 3. Renormalized diquark condensate (7) as a functionof chemical potential. Dashed red line represents the fit by(22), for detailed discussion see the paragraph after eq. (26). mass ma are completely equivalent at zero baryon densityand their multiplicative renormalization factors coincide.As the consequence, the multiplicative renormalizationof the diquark condensate (7) can be taken to be equalto the multiplicative renormalization factor of the chiralcondensate (6).To calculate the renormalized chiral and diquark con-densates one needs to know the value of the constant F . To find this constant we fit our lattice results forthe quark number density by the ChPT formula (24)in the region aµ ∈ (0 , .
12) with the pion mass m π =738 ±
13 MeV ( am π = 0 . ± . χ / dof ∼ F = 60 . ± . aF = 0 . ± . (cid:104) δ ( (cid:126)E a ) (cid:105) T = 12 N t g (cid:18) (cid:104) U tP (cid:105) µ =0 − (cid:104) U tP (cid:105) µ (cid:19) , (8) (cid:104) δ ( (cid:126)H a ) (cid:105) T = 12 N t g (cid:18) (cid:104) U sP (cid:105) µ =0 − (cid:104) U sP (cid:105) µ (cid:19) , (9)where U sP , U tP are the spatial and temporal plaquetes. Itis clear that these observables show how the chromoelec-tric and chromomagnetic fields are affected by the baryondensity.To study the topological properties of the dense two-color QCD we are going to calculate the topological sus-ceptibility. The details for these measurements mainlyfollow [43]. To get the topological charge on each config-uration we use the Gradient Flow technique [44, 45]. Onthe smoothened configurations we measure: Q L = − π (cid:88) x ± (cid:88) µνρσ = ± ˜ (cid:15) µνρσ Tr U µν ( x ) U ρσ ( x ) , (10)where U µν ( x ) is the plaquette at the point x in directions µ and ν . The final estimator for the topological charge Q is given by: Q = round ( αQ L ) , (11)where round gives the closest integer to its argument andthe factor α is chosen in such a way that it minimizes (cid:104) ( αQ L − round ( αQ L )) (cid:105) . (12)By doing this we rescale the topological charge Q L sothat its peaks become closer to integer values and thenround the estimation to this integer, thus reducing latticeartifacts. The topological susceptibility is then given by a χ top = (cid:104) Q (cid:105) N s N τ . (13)Finally we are going to study the equation of state(EoS) for the dense two-color QCD. The pressure p can becalculated if one takes the integral of the quark numberdensity p ( µ ) = (cid:90) µ dξn q ( ξ ) + p (0) . (14)In this paper we work at zero temperature, where it isreasonable to take p (0) = 0. In the calculation, thebaryon density was interpolated by cubic splines. Thepressure was then obtained with numerical integration ofthe interpolated baryon density. .
15 0 .
20 0 .
25 0 .
30 0 .
35 0 . µa . . . . . . . . . n q / n
800 1000 1200 1400 1600 1800 µ, MeV
Figure 4. The ratio n q /n as a function of chemical potential,where n = 4 µ / π is the quark number density for freerelativistic quarks. Next let us consider the trace anomaly I ( µ ) = (cid:104) T µµ (cid:105) = (cid:15) ( µ ) − p ( µ ) which can be written as [46, 47] I ( µ ) T = I G ( µ ) T + I F ( µ ) T (15) I G ( µ ) T = N t β ( g ) (cid:2) (cid:104) S G (cid:105) µ − (cid:104) S G (cid:105) (cid:3) , (16) I F ( µ ) T = − N t γ ( g ) ma (cid:2) (cid:104) ¯ qq (cid:105) µ − (cid:104) ¯ qq (cid:105) (cid:3) , (17)where I G ( µ ), I F ( µ ) denote the gluon and fermion con-tributions to the anomaly, (cid:104) S G (cid:105) is the averaged value ofthe tree level improved Symanzik gauge action, and β ( g )and γ ( g ) are β ( g ) = 4 dg − d log a , (18) d log( ma ) d log a = γ ( g ) . (19)Here two comments are in order.First, in Eq. (15) we subtracted the trace anomaly atzero temperature and density in order to get rid of theadditive divergence. Second, all simulations are carriedout on the lattice 32 . Although one cannot exclude fi-nite temperature or finite N t effects, the temperature inour simulations is close to zero. When we divide someobservable by T in the corresponding power (similar toformulae (8–9) or (15–17) ), instead of T we imply inversetemporal size of our lattice 1 /aN t (cid:39)
140 MeV.Having calculated the pressure and the trace anomalywe can calculate the energy density (cid:15) and the entropydensity s (cid:15) = I + 3 p (20) s = (cid:15) + p − µn q T (21) .
15 0 .
20 0 .
25 0 .
30 0 .
35 0 . µa . . . . . . . m π µ Σ qq
800 1000 1200 1400 1600 1800 µ, MeV
Figure 5. The ratio m π Σ qq /µ as a function of chemical po-tential, where Σ qq is defined in (7). III. FERMIONIC OBSERVABLES
The phase diagram of dense two-color QCD was stud-ied previously within the ChPT in papers [22, 23]. Thisphase diagram at zero λ can be described as follows: forsufficiently small chemical potential the chiral symmetryis broken and the chiral condensate takes nonzero val-ues, while the diquark condensate and baryon density arezero; at µ c = m π / µ > µ c the chi-ral condensate and the baryon density become nonzero.In the lattice formulation due to the finite pion massthe chiral condensate is nonzero after the transition, butit decreases with chemical potential. Moreover, nonzerovalues of the λ parameter change the second order phasetransition to a crossover. At the leading order approxi-mation of the ChPT, the dependence of the diquark con-densate, chiral condensate and quark number density onthe chemical potential can be described by the followingformulae [23] (cid:104) qq (cid:105) = 2 N f G sin α , (22) (cid:104) ¯ qq (cid:105) = 2 N f G cos α , (23) n = 8 N f F µ sin α , (24)where the α angle can be extracted from the equation µ sin α cos α = µ c (cid:18) sin α − λm cos α (cid:19) , (25)and the constant G = 12 N f (cid:112) (cid:104) ¯ qq (cid:105) + (cid:104) qq (cid:105) . (26)Our lattice results for the renormalized diquark conden-sate Σ qq , the chiral condensate Σ ¯ qq and the quark numberdensity are presented in Fig. 3, 2, 1 correspondingly.To proceed we fit simultaneously our lattice data forthe diquark condensate, the chiral condensate and thequark number density in the region aµ ∈ (0 , .
12) bymodified formulae (22)-(24). The modification consistsin addition of one constant c to the quark condensateand the constant c to the chiral condensate, these con-stants are aimed at account of the additive divergenceswhich are contained in the lattice results. Thus in thefitting procedure we have three parameters: F, µ c , c , c .In Figs. 3, 2, 1 we present the results of this fit. Fromthese figures it may be seen, that the fit quality is good( χ / dof ∼
1) and F = 63 ± µ c = 403 ±
30 MeV( aF = 0 . ± . aµ c = 0 . ± . µ c obtained in the fitting proce-dure within the uncertainty agrees with that calculatedfrom the pion mass: µ c = m π / ± µ <
540 MeV) the systemunder study is well described by ChPT. A similar con-clusion was also drawn in papers [29, 32, 33]. From Fig. 3 and Fig. 1 it may be observed that inthe region µ >
540 MeV ( aµ > .
12) the lattice data forthe diquark condensate and quark number density startto deviate from the leading order ChPT predictions. Tounderstand the origin of this deviation we remind, thatat the leading order approximation of the chiral pertur-bation theory one may ignore the interactions betweenhadrons, what can be done when baryons form a dilutegas, where the interactions are not important. It is clearthat the larger the baryon density the more importantthe interactions between baryons are, the larger the de-viation from the leading order ChPT. Thus the deviationof lattice data from ChPT predictions can be consideredas the transition of the system from a dilute baryon gasto dense matter phase.For sufficiently large density the wave functions of dif-ferent baryons overlap. If the density is increased further,an individual quark no longer belongs to a particularbaryon. One can expect that in this region of the chem-ical potentials, the system is similar in some propertiesto the Bardeen-Cooper-Schrieffer theory. Following [28]below we refer this region to as the BCS phase. In thisphase the relevant degrees of freedom are quarks forminga Fermi sphere, and the baryon density is given by theone of non-interacting quarks n = 4 µ / π . In otherwords, in the two-color QCD the diquarks are Cooperpairs in the BCS theory.Note that the notion “BCS phase” is not fully appro-priate to describe the system under study in the regionof sufficiently large baryon densities. This is because theresults of the BCS theory are applicable in the weak cou-pling regime, which might take place at ultrahigh den-sities only. In particular, one of the predictions of theBCS theory is Σ qq ∼ ∆( µ ) µ [48], where the ∆( µ ) isthe mass gap in the fermionic spectrum. In the weakcoupling regime ∆( µ ) ∼ µg − exp ( − π / √ g ) [13]. Butthe baryon densities reached in our studies are moder-ate and the system in this region is still strongly cou-pled [37], thus the weak coupling formula for the massgap is not applicable, while one might expect that therelation Σ qq ∼ ∆( µ ) µ survives. The factor ∼ µ in thelast formula results from the fermionic density of stateson the Fermi surface, and the factor ∆( µ ) determines thestrength the U V (1) symmetry breaking in the system.In order to find the value of the chemical potentialwhere the BCS phase is formed, in Figs. 4, 5 we plot theratios n q /n and m π Σ qq /µ , respectively. From Fig. 4it may be observed, that in the region aµ ∈ (0 . , . n/n deviates from the unity by not more than10 %. What concerns the ratio m π Σ qq /µ , it goes to aplateau in the same region, i.e. the condensation of di-quarks takes place on the surface of the Fermi sphereand the ∆( µ ) weakly depends on µ . From these obser-vations we can conclude, that in the region aµ > . . . . . . µa − h δ ( ~ H a ) i / T µ, MeV
Figure 6. The ratio (cid:104) δ ( (cid:126)H a ) (cid:105) /T , defined in (9), as a functionof chemical potential. IV. GLUONIC OBSERVABLES
This section is devoted to the study of the gluon back-ground at nonzero baryon density. To this end, in Fig. 7,6 we plot ratios (8), (9) as functions of baryon chemicalpotential.From Fig. 7 it may be observed, that the chromoelec-tric field decreases with increasing baryon density. Webelieve that this behaviour can be attributed to wellknown Debye screening of chromoelectric field in densematter. This phenomenon was also observed in the studyof Polyakov loop correlators [37] and gluon propaga-tors [38, 39] in dense matter. It is interesting to notethat in the BCS phase, chromoelectric field scales as −(cid:104) δ ( (cid:126)E a ) (cid:105) ∼ µ .Next let us consider the chromomagnetic field shownin Fig. 6. From this plot it is seen that within the un-certainty chromomagnetic field does not change as com-pared to its vacuum value up to aµ ∼ .
2. In the region aµ > . aµ > . i.e. spatial confinement playsa less important role, thus chromomagnetic field is lessscreened. Similar results were obtained in papers [38, 39].To study how nonzero baryon density influences thetopological properties of QC D, we calculated the topo-logical susceptibility χ for various values of the baryonchemical potentials under study. The result of this calcu-lation is presented in Fig. 8. Despite large uncertaintiesat a few points, it may be seen from this plot that thetopological susceptibility slowly decreases with rising ofthe chemical potential. This result is in disagreementwith one recent study [33], but it agrees with the resultsof the papers [49, 50]. We believe that our findings on . . . . . µa − h δ ( ~ E a ) i / T µ, MeV
Figure 7. The ratio −(cid:104) δ ( (cid:126)E a ) (cid:105) /T , defined in (8), as a func-tion of chemical potential. The minus sign is taken since thechromoelectric field decreases in dense matter as compared tothe vacuum value. the topological susceptibility are supported by the otherresults of this paper. In particular, in this section we ob-served, that chromoelectric fields are screened in densematter. For this reason one can expect that topologicalfluctuations are suppressed by dense matter as comparedto vacuum. V. EQUATION OF STATE OF DENSE QC D To study the EoS of dense two-color matter in thispaper we are going to use equations (14)-(21). For the . . . . . µa . . . . . . . χ × , G e V µ, MeV
Figure 8. Topological susceptibility in energy units, scaledby 10 for the better visual presentation, as a function of thechemical potential. .
15 0 .
20 0 .
25 0 .
30 0 .
35 0 .
40 0 . µa . . . . . .
18 1 /π (cid:15)/µ p/µ
800 1000 1200 1400 1600 1800 2000 µ, MeV
Figure 9. The energy density and the pressure divided by µ as a function of chemical potential (blue circles are slightlyshifted for the better visibility). The dashed line correspondsto the (cid:15) and 3 p of a free relativistic quark gas (cid:15) = 3 p = µ /π . functions β ( g ) the two-loop perturbative expression isused, which is independent on regularization: β ( g ) = − g (cid:18) g π + 29 g π (cid:19) . (27)For the function γ ( g ) in the calculation we use the pertur-bative one-loop expression (28), which similarly to (27)does not depend on regularization γ ( g ) = 1 + 932 π g . (28)Notice, that the use of the one-loop expression might leadto systematic uncertainty. However, we believe that weare close to the continuum limit and the expression (28) isa good approximation for the actual γ ( g ) function. Thisstatement is supported by the findings in [51] devotedto the EoS of SU (3) QCD. In this paper it was found,that for sufficiently small g the γ ( g ) is well described byone-loop formula similar to (28). The fact that we workclose to the continuum limit can be seen from the follow-ing observation: if instead of the two-loop β ( g ) functionone uses one-loop expression the results for the β ( g ) willchange by 10 %. Moreover, the γ ( g ) function enters tothe fermion contribution to the anomaly. From what fol-lows it will be clear, that the fermion contribution to theanomaly is quite small and it modifies the EoS withinthe uncertainty of the calculation.In the Fig. 10 we plot the gluon I G and fermion I F con-tributions to the anomaly and the pressure p as a functionof the chemical potential. In order to plot these observ-ables in one figure we have rescaled them. The energydensity (20) is the sum of the I G , I F and 3 p , thus usingFig. 10 one can explore the role of these three contribu-tions in the EoS. From the Fig. 10 it may be seen that the .
15 0 .
20 0 .
25 0 .
30 0 .
35 0 .
40 0 . µa I G ( µ ) /T I F ( µ ) /T p/ T
800 1000 1200 1400 1600 1800 2000 µ, MeV
Figure 10. The gluon I G and fermion I F contributions to theanomaly, defined in (16) and (17) respectively, and pressure p as functions of the chemical potential. In order to plot theseobservables in one figure we rescaled them. smallest contribution is the fermion part of the anomaly, I F . The next (by the size) is the gluon contribution tothe anomaly, I G . Unfortunately, the uncertainty of thecalculation of this observable is quite large, and the I G develops nonzero values only in the region aµ ≥ .
4. No-tice that for all the values of µ under study the fermioniccontribution I F is much smaller than the uncertainty ofthe calculation of the gluonic contribution I G . The term3 p provides the largest contribution to the energy den-sity (20). In the region, which can be well describedby ChPT, aµ < .
12, the 3 p term is compatible to theuncertainty in I G . In the following region of the phasediagram, where the system becomes dense, aµ > . p term is larger than the gluonic contribution I G ,but the uncertainty in the energy density remains quitelarge. Finally, in the BCS phase the uncertainty in theenergy density (cid:15) becomes small.In Fig. 9 we present the energy density and the pres-sure, divided by µ , as functions of the chemical poten-tial. The dashed line corresponds to the EoS of free rel-ativistic quark gas (cid:15) = 3 p = µ /π . It may be observedfrom the Fig. 9, that in the BCS phase the EoS is welldescribed by the EoS of free relativistic quarks. In oursimulations the quark mass is quite large leading to largepion mass, but nonzero quark mass does not play an im-portant role in the BCS phase.Now let us focus on the entropy density (21). Ac-cording to the third law of thermodynamics the entropyapproaches to some constant value as temperature ap-proaches to zero. If the ground state is not degenerate,the entropy is zero. In our simulations due to nonzerovalue of the λ -parameter there is no degeneracy in thesystem under study, thus one can expect that the en- μ m π /2 hadronic phase deconfinementdense quark matterBEC BCS
900 MeV
Figure 11. Schematic phase diagram of dense two-color QCDat low temperatures. tropy is zero. Our results confirm that s = 0 within theuncertainty of the calculation for all values of the chem-ical potential under consideration.Previously the EoS of dense QC D was studied in thepapers [28, 34], where lattice simulation was carried outwith dynamical Wilson fermions. It is rather difficult tocompare our results and the results obtained there dueto large uncertainties of the calculation at small valuesof chemical potential. However, in the BCS phase theEoS from [28, 34] is well described by the EoS of freerelativistic quarks, which agrees with the results of thepresent paper.
VI. DISCUSSION AND CONCLUSION
In this paper, we carried out lattice study of the phasediagram of dense two-color QCD with N f = 2 quarks andthermodynamic properties of this system. This studywas conducted at low temperature and for the baryonchemical potential in the region µ ∈ (0 , µ < m π /
2) the system isin the hadronic phase, and the chiral symmetry is bro-ken; at µ = m π / µ ≥ m π /
2, but not at too large chemical poten-tial, is similar to a dilute baryon gas. Lattice results forthe baryon density, diquark and chiral condensates arewell described by ChPT up to µ <
540 MeV.Increasing the baryon density further, we proceed todense matter, where the interactions between baryonscannot be fully accounted within perturbation theory. This transition manifests itself in terms of the deviationof different observables from the ChPT predictions. Inparticular, in this paper the deviation is well pronouncedin the diquark condensate and the baryon density.At sufficiently large baryon density ( µ ∼
900 MeV, aµ ∼ .
20) some observables of the system under studycan be described using Bardeen-Cooper-Schrieffer theory(BCS phase). In particular, the baryon density is well de-scribed by the density of non-interacting fermions whichoccupy a Fermi sphere of radius r F = µ . Moreover, thediquark condensate, which plays the role of a condensateof Cooper pairs, is proportional to the Fermi surface.In the region aµ < . µ ∼
900 MeV( aµ ∼ .
2) we observe confinement/deconfinement tran-sition in dense two-color QCD [36]. This transition man-ifests itself in a rise of the Polyakov loop and vanishingof the string tension. It was also found that, after decon-finement is achieved, spatial string tension σ s decreasesmonotonically and ends up vanishing at µ q ≥ aµ ≥ . aµ ∼ .
2, then in the region aµ > . aµ > . i.e. spatial confinement plays less importantrole. As the result chromomagnetic field is less screened.Similar results were obtained in papers [38, 39].To study how nonzero baryon density influences thetopological properties of QC D we have calculated thetopological susceptibility χ . The topological susceptibil-ity slowly decreases with rising of the chemical potential.We believe, that this decrease of χ with increasing chem-ical potential is related to the screening of chromoelectricfields in dense matter.Finally in this paper the equation of state of denseQC D was studied. In particular, we calculated the pres-sure, the trace anomaly and the energy density. Althoughit is possible to calculate the pressure with rather goodaccuracy, the uncertainty in the trace anomaly is ratherlarge. As a result, good accuracy in the energy densitycan be achieved only in the BCS phase. It is interest-ing to note, that in the BCS phase the equation of stateis well described by the corresponding formulae for freerelativistic fermions (cid:15) (cid:39) p, p (cid:39) µ / π . The entropydensity remains zero within the uncertainty of the calcu-lation for all values of the chemical potential as it shouldbe in the λ (cid:54) = 0 case.At the end of this paper we are going to discuss thedependence of results on the diquark source λ . Withinthe accuracy of the calculation we don’t see any λ -dependence of the gluonic observables. As for thefermionic observables in the region of the phase dia-gram, which is well described by the ChPT, there is λ -dependence. In particular, both the quark number den-sity and the diquark condensate are zero for λ = 0 inthe region µ < m π /
2, and finite value of the λ -parameterleads to nonzero values of these observables. Close to thephase transition µ ∼ m π / λ -dependenceof the fermionic observables. However, when one movesfrom the region of the phase transition µ ∼ m π / λ -dependence becomes weak andthe fermionic observables at λ/m ∼ . λ →
0. This result is in agreementwith the previous findings in [29, 32]. The λ -dependenceof the EoS is related to the mentioned λ -dependence ofthe quark number density, because the chiral condensatedoes not depend on it. In the region µ < m π / µ ∼ m π / λ is weak, thus the λ -dependence of the EoS is weak. ACKNOWLEDGMENTS
We would like to thank Prof. Naoki Yamamoto foruseful discussions. The work of V. V. B., which con-sisted of study of the equation of state and physicalinterpretation of lattice data, was supported by grantfrom the Russian Science Foundation (project num-ber 16-12-10059). N. Yu. A., A. Yu. K. and A. A. N. ac-knowledge the support of the Russian Foundation forBasic Research (project no. 18-32-20172 mol a ved).A. A. N. is also grateful for the support from STFCgrant ST/P00055X/1. This work has been carried outusing computing resources of the federal collective us-age center Complex for Simulation and Data Process-ing for Mega-science Facilities at NRC “Kurchatov Insti-tute”, http://ckp.nrcki.ru/ . In addition, the authorsused the cluster of the Institute for Theoretical and Ex-perimental Physics and the supercomputer of Joint In-stitute for Nuclear Research “Govorun”. [1] S. Muroya, A. Nakamura, C. Nonaka, and T. Takaishi,Prog. Theor. Phys. , 615 (2003), arXiv:hep-lat/0306031 [hep-lat].[2] C. Ratti, PoS
LATTICE2018 , 004 (2019).[3] M. G. Alford, K. Rajagopal, and F. Wilczek, Nucl. Phys.
B537 , 443 (1999), arXiv:hep-ph/9804403 [hep-ph].[4] M. G. Alford, J. Berges, and K. Rajagopal, Nucl. Phys.B , 219 (1999), arXiv:hep-ph/9903502.[5] M. Buballa, Phys. Rept. , 205 (2005), arXiv:hep-ph/0402234 [hep-ph].[6] T. Khunjua, K. Klimenko, and R. Zhokhov, Phys. Rev.D , 054036 (2018), arXiv:1710.09706 [hep-ph].[7] T. Khunjua, K. Klimenko, and R. Zhokhov, Phys. Rev.D , 054030 (2018), arXiv:1804.01014 [hep-ph].[8] T. G. Khunjua, K. G. Klimenko, and R. N. Zhokhov,Particles , 62 (2020).[9] Z. Khaidukov and Y. A. Simonov, Phys. Rev. D ,076009 (2019), arXiv:1906.08677 [hep-ph].[10] C. S. Fischer, Prog. Part. Nucl. Phys. , 1 (2019),arXiv:1810.12938 [hep-ph].[11] W.-j. Fu, J. M. Pawlowski, and F. Rennecke, Phys. Rev.D , 054032 (2020), arXiv:1909.02991 [hep-ph].[12] F. Gao and J. M. Pawlowski, (2020), arXiv:2002.07500[hep-ph].[13] D. T. Son, Phys. Rev. D59 , 094019 (1999), arXiv:hep-ph/9812287 [hep-ph].[14] L. McLerran and R. D. Pisarski, Nucl. Phys.
A796 , 83(2007), arXiv:0706.2191 [hep-ph].[15] T. Gorda and P. Romatschke, Phys. Rev. D , 014019(2015), arXiv:1412.6712 [hep-ph].[16] T. Gorda, A. Kurkela, P. Romatschke, M. Sppi, andA. Vuorinen, Phys. Rev. Lett. , 202701 (2018), arXiv:1807.04120 [hep-ph].[17] D. T. Son and M. A. Stephanov, Phys. Rev. Lett. ,592 (2001), arXiv:hep-ph/0005225 [hep-ph].[18] J. Kogut and D. Sinclair, Phys. Rev. D , 034505 (2002),arXiv:hep-lat/0202028.[19] J. Kogut and D. Sinclair, Phys. Rev. D , 014508 (2002),arXiv:hep-lat/0201017.[20] B. B. Brandt, G. Endrodi, and S. Schmalzbauer, Phys.Rev. D97 , 054514 (2018), arXiv:1712.08190 [hep-lat].[21] B. B. Brandt, G. Endrodi, and S. Schmalzbauer, PoS
Confinement2018 , 260 (2018), arXiv:1811.06004 [hep-lat].[22] J. B. Kogut, M. A. Stephanov, and D. Toublan, Phys.Lett.
B464 , 183 (1999), arXiv:hep-ph/9906346 [hep-ph].[23] J. B. Kogut, M. A. Stephanov, D. Toublan, J. J. M.Verbaarschot, and A. Zhitnitsky, Nucl. Phys.
B582 , 477(2000), arXiv:hep-ph/0001171 [hep-ph].[24] S. Hands, J. B. Kogut, M.-P. Lombardo, and S. E.Morrison, Nucl. Phys.
B558 , 327 (1999), arXiv:hep-lat/9902034 [hep-lat].[25] J. B. Kogut, D. Toublan, and D. K. Sinclair, Phys. Lett.
B514 , 77 (2001), arXiv:hep-lat/0104010 [hep-lat].[26] J. B. Kogut, D. Toublan, and D. K. Sinclair, Nucl. Phys.
B642 , 181 (2002), arXiv:hep-lat/0205019 [hep-lat].[27] S. Muroya, A. Nakamura, and C. Nonaka, Phys. Lett.
B551 , 305 (2003), arXiv:hep-lat/0211010 [hep-lat].[28] S. Cotter, P. Giudice, S. Hands, and J.-I. Skullerud,Phys. Rev.
D87 , 034507 (2013), arXiv:1210.4496 [hep-lat].[29] V. V. Braguta, E. M. Ilgenfritz, A. Yu. Kotov, A. V.Molochkov, and A. A. Nikolaev, Phys. Rev.
D94 , 114510(2016), arXiv:1605.04090 [hep-lat]. [30] L. Holicki, J. Wilhelm, D. Smith, B. Wellegehausen, andL. von Smekal, Proceedings, 34th International Sympo-sium on Lattice Field Theory (Lattice 2016): Southamp-ton, UK, July 24-30, 2016 , PoS
LATTICE2016 , 052(2017), arXiv:1701.04664 [hep-lat].[31] T. Boz, O. Hajizadeh, A. Maas, and J.-I. Skullerud,Phys. Rev. D , 074514 (2019), arXiv:1812.08517 [hep-lat].[32] J. Wilhelm, L. Holicki, D. Smith, B. Wellegehausen,and L. von Smekal, Phys. Rev. D , 114507 (2019),arXiv:1910.04495 [hep-lat].[33] K. Iida, E. Itou, and T.-G. Lee, JHEP , 181 (2020),arXiv:1910.07872 [hep-lat].[34] T. Boz, P. Giudice, S. Hands, and J.-I. Skullerud, Phys.Rev. D , 074506 (2020), arXiv:1912.10975 [hep-lat].[35] P. Buividovich, D. Smith, and L. von Smekal, (2020),arXiv:2007.05639 [hep-lat].[36] V. G. Bornyakov, V. V. Braguta, E. M. Ilgenfritz, A. Yu.Kotov, A. V. Molochkov, and A. A. Nikolaev, JHEP ,161 (2018), arXiv:1711.01869 [hep-lat].[37] N. Astrakhantsev, V. Bornyakov, V. Braguta, E.-M. Il-genfritz, A. Kotov, A. Nikolaev, and A. Rothkopf, JHEP , 171 (2019), arXiv:1808.06466 [hep-lat].[38] V. Bornyakov, A. Kotov, A. Nikolaev, and R. Rogalyov,Particles , 308 (2020), arXiv:1912.08529 [hep-lat]. [39] V. Bornyakov, V. Braguta, A. Nikolaev, and R. Roga-lyov, (2020), arXiv:2003.00232 [hep-lat].[40] P. Weisz, Nucl. Phys. B212 , 1 (1983).[41] G. Curci, P. Menotti, and G. Paffuti, Phys. Lett. ,205 (1983), [Erratum: Phys. Lett.135B,516(1984)].[42] A. Bazavov et al. , Phys. Rev.
D85 , 054503 (2012),arXiv:1111.1710 [hep-lat].[43] C. Bonati and M. D’Elia, Phys. Rev. D , 105005(2014), arXiv:1401.2441 [hep-lat].[44] M. Luscher, Commun. Math. Phys. , 899 (2010),arXiv:0907.5491 [hep-lat].[45] M. Lscher, JHEP , 071 (2010), [Erratum: JHEP 03,092 (2014)], arXiv:1006.4518 [hep-lat].[46] S. Borsanyi, G. Endrodi, Z. Fodor, A. Jakovac, S. D.Katz, S. Krieg, C. Ratti, and K. K. Szabo, JHEP ,077 (2010), arXiv:1007.2580 [hep-lat].[47] S. Borsanyi, G. Endrodi, Z. Fodor, S. Katz, S. Krieg,C. Ratti, and K. Szabo, JHEP , 053 (2012),arXiv:1204.6710 [hep-lat].[48] T. Kanazawa, T. Wettig, and N. Yamamoto, JHEP ,003 (2009), arXiv:0906.3579 [hep-ph].[49] B. Alles, M. D’Elia, and M. Lombardo, Nucl. Phys. B , 124 (2006), arXiv:hep-lat/0602022.[50] B. Alles, M. D’Elia, and M. Lombardo, Nucl. Phys. BProc. Suppl. , 225 (2007), arXiv:hep-lat/0609017.[51] M. Cheng et al. , Phys. Rev. D77