Probing QCD critical point and induced gravitational wave by black hole physics

Locating the critical endpoint of QCD and the region of a first-order phase transition at finite baryon chemical potential is an active research area for QCD matter. We provide a gravitational dual description of QCD matter at finite baryon chemical potential and finite temperature using the non-perturbative approach from gauge/gravity duality. After fixing all model parameters using state-of-the-art lattice QCD data at zero chemical potential, the predicted equations of state and QCD trace anomaly relation are in quantitative agreement with the latest lattice results. We then give the exact location of the critical endpoint as well as the first-order transition line, which is within the coverage of many upcoming experimental measurements. Moreover, using the data from our model at finite baryon chemical potential, we calculate the spectrum of the stochastic gravitational wave background associated with the first-order QCD transition in the early universe, which could be observable via pulsar timing in the future.

Introduction-As one of the most interesting and fundamental challenges of high energy physics, the phase diagram of QCD has been intensively studied. It involves the behavior of strongly interacting matter under extreme conditions [1][2][3], ranging from cosmology and astrophysics to heavy-iron collisions. Due to the strong interaction in the non-perturbative regime, it has not been possible to obtain the full QCD phase diagram directly from QCD in terms of the temperature T and the baryon chemical potential µ B or the baryon number density n B . While lattice QCD can provide reliable first-principles information at zero density [4][5][6], it fails at finite density due to the famous sign problem [2]. Note, however, that the lattice data can be extrapolated to finite µ B via different systematical schemes [7][8][9], but this is only reliable for small µ B . Meanwhile, many effective low energy models have been proposed to study the QCD phase diagram under certain conditions, such as [10][11][12], for which it is difficult to match the lattice QCD data quantitatively. Nevertheless, the QCD matter is now believed to be in a hadronic phase of color-neutral bound states at low T and small µ B , while it is deconfined at high T and large µ B , known as the quark-gluon plasma. These two phases are separated by a transition at small µ B and are expected to change into a first-order transition for higher µ B . The critical point between them is the QCD critical endpoint (CEP), which has been under active investigation by experiments and theoretical calculations.
On the other hand, the gauge/gravity duality [13][14][15][16] provides a powerful non-perturbative approach to solving the strongly coupled non-Abelian gauge theories by mapping to a weakly coupled gravitational system with one higher dimension. In particular, it provides a convenient way to incorporate real-time dynamics and transport properties at finite temperatures and densities. Previous studies (e.g. [17][18][19][20][21]) have provided a strong indication that holography can make quantitative description of the properties of QCD in the non-perturbative regime. The construction of the QCD phase diagram in such a holographic approach was initiated by [22,23], where the Einstein-Maxwell-dilaton (EMD) theory was used to mimic properties in the QCD phase diagram. Since then, several attempts have been made toward the direction of the phase diagram in the T − µ B plane (see e.g. [24][25][26][27][28][29][30][31][32]). In the spirit of effective field theory, the model parameters of the bulk gravitational theory should be fixed by matching with lattice QCD results. Therefore, it is crucial to use up-to-date lattice simulation for reliable prediction at finite µ B .
In this Letter, we construct a holographic QCD (hQCD) model in which all parameters are fixed using state-of-the-art lattice QCD data [6] at µ B = 0, generated by highly improved stagger fermion action. Moreover, all thermodynamic quantities are computed directly from the holographic renormalization and the so-called thermodynamic consistency relations [33]. Our prediction for the thermodynamic observables at finite µ B is in quantitative agreement with the latest lattice results that are available for µ B /T < 3.5 [9]. We also calculate both the chiral and gluon condensates. Remarkably, we show for the first time in holography that the gluon condensate agrees with the lattice simulation regarding the QCD conformal anomaly. We then manage to make precise predictions for the QCD phase diagram at finite arXiv:2201.02004v4 [hep-th] 18 Feb 2023 µ B , in particular a first-order line and a CEP located at T C = 105MeV and µ C = 555MeV. Interestingly, the location of our CEP can be checked in the near future by many important upcoming facilities. Moreover, the strong first-order phase transition (SFOPT) in the early universe is an important source of the stochastic gravitational wave (GW) background (see, e.g. [38][39][40][41] and references therein). While various scenarios of new physics beyond the Standard Model of particle physics have been considered to engineer an SFOPT in the literature, our present model provides a scenario for phase transition GWs within the Standard Model.
Holographic model -We begin with the fivedimensional EMD theory with a minimal set of fields for capturing the essential dynamics. Here κ 2 N is the effective Newton constant. In addition to the metric g µν that characterizes the geometry, the real scalar φ (known as dilaton) encodes the running of the gauge coupling, and the Maxwell field A µ accounts for a finite baryon density. Z(φ) and V (φ) are two phenomenological terms that will be fixed by matching to the lattice QCD at µ B = 0.
The hairy black hole reads where dx 2 3 = dx 2 + dy 2 + dz 2 and r is the holographic radial coordinate with the asymptotical anti-de Sitter (AdS) boundary at r → ∞. Denoting the event horizon as r h at which f vanishes, the temperature and entropy density can be obtained as while the baryon chemical potential µ B and density n B can be obtained from A t at the AdS boundary, we read the energy density and pressure P directly by the dual stress-energy tensor via the holographic renormalization [42,43]. We refer to the Supplemental Material (SM) for more details [44]. Then, the equation of state (EOS) and transport properties can be determined precisely. The form of V and Z is partially motivated by the one in [17,22,23], although these models were not able to simultaneously fit the lattice data for equilibrium and near-equilibrium features quantitatively.
By global fitting the state-of-the-art lattice data [6,9] with (2+1)-flavors at zero net-baryon density (see Fig. 1 [6,9], while solid red curves are from our hQCD model. Nτ is the temporal extent in lattice QCD. Following lattice QCD notations, the condensates of the light (l = u, d) quarks and the s quark are normalized by there values at zero temperature. The gluon condensate T is subtracted by its vacuum value, where β(g) is the β-function and g is QCD gauge coupling. the hQCD model can be fixed to be with c 1 = 0.7100, c 2 = 0.0037, c 3 = 1.935, c 4 = 0.085, c 5 = 30. Here φ s = rφ| r→∞ is the source term of the scalar field φ, which essentially breaks the conformal symmetry and plays the role of the energy scale.
The fitting results are presented in the upper panel of Fig. 1. In addition to the EOS, the specific heat C V = (d /dT ) µ B , and the second-order baryon susceptibility χ B 2 = (dn B /dµ B ) T /T 2 at zero density, which are important quantities characterizing QCD transition, agree with the lattice data quantitatively.
Given the local bulk description of QCD (4), we also calculate the gluon condensate and the condensates for u, d, s quarks in the spirit of the KKSS model [45] by introducing appropriate probe actions. The temperature dependence of these condensates is in quantitative agreement with the lattice simulations [6], see the lower panel of Fig. 1. Remarkably, the trace anomaly ( −3P ) is found to be consistent with the summation of the quark and gluon condensates. It is the first time that this relation is realized in a holographic setup (see SM for technical details).
So far, our holographic model is completely fixed, and there are no free parameters. In Fig. 2, we compare the holographic results with the latest lattice QCD [9] which combines the Taylor-expanded approach and the analytical continuation approach. The holographic predictions are in quantitative agreement with the lattice results available for small chemical potentials, which strongly supports our hQCD model.
QCD phase diagram-Having completely fixed the holographic model, we can now construct the phase diagram in the T − µ B plane by numerically solving a series of black holes. We find that, as µ B increases, the crossover on the T -axis is sharpened into a first-order line at the CEP (see Fig. 3). Since the transition from µ B = 0 up to the CEP is a smooth crossover, there is no unique way to determine the transition temperature in the literature. Some suitable probes characterizing the drastic change of degrees of freedom between the quarkgluon plasma (QGP) and the hadron resonances gas are the minimum of the speed of sound c s and the maximum increasing point of χ B 2 . The transition lines of the two probes are shown in Fig. 3. The transition temperatures of the probes at zero µ B compare well with the predictions of lattice QCD for the same up to the quark-hadron transition region around 140−160MeV [46][47][48]. Although they do not coincide quantitatively in the crossover region, they show similar behavior and are in the same order of magnitude. Moreover, they come together at the critical point µ C .  [51][52][53][54][55][56]. NJL: PNJL is from [57] and NJL1 from [58]. FRG: F1-F2 are from [59,60]. The coverage of RHIC [61], the STAR fixed target program (FXT), and future (FAIR, JPARC-HI, and NICA) experimental facilities [65] is also indicated at the top of the figure.
For a baryon chemical potential above µ C , the QCD transition becomes a first-order transition for which the free energy can uniquely determine the transition point from the holographic calculation (see the solid black line in Fig. 3). The first order line decreases monotonically with µ B . Although we cannot see when the line will terminate, our numerics suggest that the transition temperature becomes significantly small, approximately larger than µ B ∼ 1050MeV, beyond which no stable black holes have been found. In this region, color superconductivity or related phenomena might be set in.
The CEP predicted by our hQCD model is located at T C = 105MeV, µ C = 555MeV (the red dot in Fig. 3). Therefore, there is no first-order phase transition in the region with µ B /T < µ C /T C ≈ 5.3, which is consistent with the thermodynamics of Fig. 2. We also show the location of CEP from different models and note that the current best lattice estimate suggests that CEP is likely above µ B ∼ 300MeV [49,50]. The variation in priorities is considerable, but most of them are close to the transition line we found. Moreover, compared to other predictions, our CEP is close to those obtained by the Schwinger-Dyson equation [55,56] or the functional renormalization group [59] respectively. Nevertheless, all our predictions are at least qualitatively consistent with the consensus expectations for the QCD phase diagram. Interestingly, the predicted location of CEP is within the coverage of the STAR fixed target program (FXT) and future (FAIR [62,63], JPARC-HI and NICA [64]) experimental facilities [65]. Therefore, our prediction can be verified in the near future.
GWs from QCD phase transition-An important feature associated with the first order QCD transition in the early universe is the generation of a stochastic background of GWs, for which bubbles of the broken phase nucleate and expand in the presence of a plasma background of false vacuum. The GWs can be generated in three processes: bubble collisions, sound waves, and magnetohydrodynamic turbulence. Some preliminary studies of the GW spectra from holographic QCD have been summarized in [66]. It has been recognized that higher-order corrections prevent true runaway transitions from occurring [67], thus the bubble wall terminal velocity v w < 1. For these non-runway bubbles, the GWs will be dominated by sound waves for which the energy spectrum reads [68] with the peak frequency Hz .

(6)
The strength of the GW spectrum depends on many parameters that can be obtained from our holographic model. We discuss each of these quantities in turn.
T n is approximately by the critical temperature of the first-order transition, see the solid black line in Fig. 3. Beyond the bag model approximation, the phase transi-tion strength α is given by [39] between the false (+) and true (−) vacuums, where θ = − 3P is the trace anomaly, and w = + P is the enthalpy. Note that the numerator of (7) is the latent heat for the first-order QCD phase transition at T n . The effective number of degrees of freedom g n = 45s + /(2π 2 T 3 n ). We fix v w to the typical values v w = 0.95 for which the kinetic energy efficiency coefficient κ = α 0.73+0.083 √ α+α . The only free parameter is the inverse time duration of the phase transition β/H n . Measurements of primordial element abundances and the anisotropy spectrum of the Cosmic Microwave Background (CMB) yield a strong constraint on the baryon sector, characterized by the baryon asymmetry η B ≡ n B /s. The observation value η ob B ≈ 10 −10 [69]. This quantity is conserved except during baryogenesis and a first-order phase transition. Thus, to apply to the cosmological QCD phase transition, the baryon asymmetry after the first-order transition in our model should be compatible with η ob B in the later evolution of the universe. Near the CEP, our model yields η B ∼ 10 −2 , which is eight orders of magnitude larger than the observation. Nevertheless, the value of η B for the true vacuum at the first-order phase transition line decreases significantly as µ B increases. In particular, around µ B = 1000Mev, η B reaches the same order of magnitude as the cosmological observation.
We therefore show the GW energy spectrum for µ B = 1000Mev (µ B /T ≈ 20) in Fig. 4, for which the cos-mological QCD phase transition is first order and does not violate the constraint of a tiny baryon asymmetry. This transition point is very close to the CEP by NJL model [58] (the green triangle at the bottom right of Fig. 3). To account for the uncertainty in the duration β/H n , we will scan the space with 2 ≤ β/H n ≤ 80. The energy density of GWs can reach 10 −9 around 10 −7 Hz. While out of reach of the Laser Interferometer Space Antenna (LISA), it can be detected by the International Pulsar Timing Array (IPTA) and the Square Kilometer Array (SKA). The detection from the North American Nanohertz Observatory for Gravitational Wave (NANOGrav) is possible for extreme scenarios with small values of β/H n ∼ O(1).
Moreover, the QCD first-order phase transition with µ B < 1000Mev could also be realized in the early universe by considering the "little inflation" [70,71] (see e.g. [72][73][74] for earlier work and the SM). In this scenario, the large baryon asymmetry, which can be generalized by the well-established Affleck-Dine baryogenesis [75] is diluted by a short period of inflation in which the universe is trapped in a false vacuum state of QCD. We find that the e-folding number sufficient to dilute the baryon number in the little inflation scenario is at most 7. Another scenario for a first-order QCD transition in agreement with observation is to consider a large lepton asymmetry [76,77], which we leave for future work.
Discussion-We have built a holographic model to confront the phase diagram in (2+1)-flavor QCD. The thermodynamics behaviors (EOS, transport, and condensates) are quantitatively matched with the latest lattice QCD simulation. The model has captured the characteristic QCD transition properties and offered a reliable first-order transition line and CEP in the QCD phase diagram at finite µ B . The predicted location of CEP and the GW energy spectrum of SFOPT from our hQCD model is within the detectivity of upcoming experimental facilities and, therefore, could be verified in the future. Moreover, our EOS at finite density could be crucial for supporting the experimental programs towards QCD matter, e.g. heavy-ion collisions.
Dedicating to a precise characterization of QCD matter, particularly the properties and differences of the quark-gluon plasma and hadronic phases along the first-order phase transition, is an interesting direction for further study. We have set up the preliminary hQCD model to investigate the phase transition in the QCD phase diagram quantitatively. Many other relevant physical quantities should be considered to complete the phase diagram. Moreover, the present study should be embedded in the framework of a more general and multidimensional view of the QCD phase diagram, including magnetic field, isospin, and rotation. It will be interesting to consider real-time dynamics far from equilibrium.

APPENDICES
Equations of motion-From the action (1) in the main text we can derive the equations of motion that are given as follows.
As the dual system lives in a spatial plane, we choose the Poincaré coordinates with r the radial direction in the bulk. The metric ansatz reads together with We denote the event horizon as r h at which f vanishes. Then the temperature and entropy density are given by Substituting the ansatz into (8), we obtain the following independent equations of motion.
where the prime denotes the derivative with respect to r. Both Z(φ) and V (φ) will be determined by matching to the lattice QCD at µ B = 0.
In what follows we will specify these two functions as where c 1 , c 2 , c 3 , c 4 , c 5 are free parameters of the model. These model parameters capture the physical properties of realistic QCD,e.g. EOS and baryon susceptibility. Near the AdS boundary r → ∞ where φ → 0, one has Therefore, the cosmological constant is given by Λ = −6 (the AdS radius L = 1).
To obtain the numerical solutions for f (r), η(r), φ(r) and A t (r), we should specify suitable boundary conditions at both the event horizon r h and the AdS boundary r → ∞. The smoothness of the event horizon yields the following analytic expansion in terms of (r − r h ) in the IR: After substituting (15) into the equations of motion (12), one finds four independents coefficients (r h , a h , η 0 h , φ 0 h ). On the other hand, near the AdS boundary, we obtain the following asymptotic expansion: Note that we have taken the normalization of the time coordinate at the boundary such that η(r → ∞) = 0. φ s is the source of the scalar operator of the boundary theory, which essentially breaks the conformal symmetry and plays the role of the energy scale. Before proceeding, we point out that the equations of motion (12) have two independent scaling symmetries: with λ t and λ r constants. Thus, we can first set η 0 h = 0 and r h = 1 for performing numerics. After obtaining the numerical solutions, we should use the first symmetry to satisfy the asymptotic condition η(r → ∞) = 0 and use the second one to fix the energy scale φ s . Thermodynamics-We now compute the free energy density Ω which is identified as the temperature T times the renormalized action in the Euclidean signature. Since we consider a stationary problem, the Euclidean action is related to the Minkowski one by a minus sign. Moreover, we should include the Gibbons-Hawking boundary term for a well-defined Dirichlet variational principle and a surface counterterm for removing divergence. Therefore, we have with V = dxdydz and t ∈ [0, 1/T ]. We work in the grand canonical ensemble for which the baryon chemical potential is fixed. For the model (13) we  [6,9], while solid curves from our hQCD model, including entropy density s, trace anomaly ( − 3P ), pressure P , specific heat CV , squared sound speed c 2 s and baryon susceptibility.
are considering, following the the holographic renormalization [42,43], the boundary terms take the form where h µν is the induced metric at the AdS boundary and K µν is the extrinsic curvature defined by the outward pointing normal vector to the boundary.
Employing the equations of motion (12) and the asymptotical expansion (16), we obtain The energy-momentum tensor of the dual boundary theory reads Substituting (16), we obtain with vanishing non-diagonal components. One immediately finds that P = −Ω, which is expected from thermodynamics. The trace of the energy-momentum tensor is given by It is manifest that there is a trace anomaly in the presence of the source term φ s that breaks the conformal symmetry. By integrating the second equation of (12), one has where the constant n B is nothing but the charge density that can be computed using the standard holographic dictionary. Another useful radially conserved quantity reads [79,80] which connects horizon to boundary data. Evaluating at the horizon, we obtain and therefore Q = 0 signals the extremity. On the other hand, evaluating Q at the AdS boundary, we find Therefore, we immediately obtain the expected thermodynamical relation where Ω is the free energy density (21).   [6], while blue curves are from the holographic model of [24]. The entropy density , trace anomaly ( − 3P ), and pressure P fail to fit the lattice data for the temperature above 280MeV.
After obtaining the above thermodynamic quantities, we can also compute some important transport coefficients, including the sound speed C s = (dP/d ) µ B , the specific heat C V = (d /dT ) µ B , and the second-order baryon susceptibility χ B 2 = (dn B /dµ B ) T /T 2 . These properties are compared to the state-of-the-art lattice data [6,9] with (2 + 1)-flavors at zero baryon density, from which all free parameters of our hQCD model can be fixed, see Eq. (4) and Figs. 1 and 2 in the main text. Moreover, the parameter b that appears in the boundary terms (20) is chosen to be b = −0.27341. We stress that b is vital to build a hQCD model, which is necessary to satisfy the lattice QCD simulation for which P (T = 0) = 0 at zero baryon chemical potential. Without a concrete holographic renormalization, one cannot fix b. As far as we know, this important point has been overlooked in the literature. In most studies, the thermodynamic variables were obtained only by integrating the standard first law of black hole thermodynamics in a finite temperature range. We show the fitting results for more quantities in Fig. 5. In addition to the EOS (left panel), the sound speed c s = (dP/d ) µ B , the specific heat C V = (d /dT ) µ B , and the second-order baryon susceptibility χ B 2 = (dn B /dµ B ) T /T 2 at zero density, which are three important quantities characterizing QCD transition, agree with the lattice data quantitatively, see the right panel of Fig. 5.
In the real world, the quarks are massive, the chiral symmetry is not an exact symmetry of QCD, and the center symmetry is not exact. There are no well-established order parameters to probe the phase transition. It is nontrivial to identify whether the phase transition is associated with the chiral phase transition or deconfinement phase transition. Nevertheless, the defining feature for a first-order phase transition can be identified from the free energy computed directly in our holographic model, see (21) and (29). In Fig. 6, we show the free energy density Ω and entropy density s as a function of T for different µ B . The temperature dependence of the free energy decreases smoothly for µ B < µ C . At the same time, it becomes a swallowtail for µ B > µ C , singling a first-order phase transition that has been widely used to mimic the phase transition in holographic QCD. The location of the CEP is found to be µ c = 555MeV and T c = 105MeV. The corresponding behavior of the entropy density versus T is presented in the right panel of Fig. 6. A first-order jump in the entropy density manifests when µ B > µ C .
Before ending this section, we compare the thermodynamics of the holographic model [24] to the state-of-the-art lattice QCD data [6] used in the present work. This holographic model was constructed to mimic the QCD thermodynamics at a quantitative level. As can be seen from Fig. 7, the model of [24] fails to fit the new lattice data for the temperature above 280MeV at µ B = 0. Moreover, it is much more challenging to fit n B . One can see from the left panel of Fig. 8 that the model of [24] can only quantitatively match the latest QCD results for baryon density below µ B /T = 2.0. This fitting improves significantly in our model, see the right panel. Moreover, our setup is much simpler than the one in [24], which, as far as we know, was the only effective holographic model that agreed with the old lattice QCD data [5] for EOS and transport coefficients at the quantitative level.
Quark and gluon condensates-The nonperturbative quantum fluctuations generate quark and gluon condensates, which break chiral or dilatation symmetries of QCD. While there have been intensive studies in the literature, a complete and non-perturbative understanding of realistic QCD is still missing, and many features of these condensates are not yet well understood and established. Our holographic model has been shown to describe QCD thermodynamics quite well. A natural question is how to incorporate the flavor dynamics and gluon condensate in our framework. Following the spirit of the KKSS model [45], we introduce a holographic probe action to evaluate the chiral and gluon condensates.
We consider the (2+1)-flavor QCD for which m u = m d < m s with m u , m d and m s the mass for u, d and s quarks, respectively. Treating the flavor part as a probe, we introduce the effective form of the action in the background (9), where q = l denotes the light quarks (u and d) and q = s is for the s quark. The bulk field X q is dual to the boundaryψ q ψ q chiral operator with the scaling dimension ∆ c = 3. The effective coupling Z q (φ) and the potential V (X q ) will be fixed later by matching the lattice QCD data. Then, we need to solve the equation of motion derived from (30).
Since the dynamics of u, d and s quarks in 2+1 flavored QCD is quite different, we introduce Z l (φ) and Z s (φ) to distinguish the dynamics of light and heavy quarks. We choose Z q (q = l, s) and V as follows with a 0 , a 1 , a l , a s constants. Note that the mass of the bulk field X q is determined by the scaling dimension of the dual operatorψ q ψ q . Therefore, the UV expansion could be obtained as where m l = 4.5MeV and m s = 90MeV are the light and heavy quark masses, respectively [6]. According to the holographic dictionary, the chiral condensates read ψ ψ q,T = δS ren By requiring the regularity of the solution at the event horizon, one can obtain the chiral condensates for light and heavy quarks in the function of T . The lattice simulation [81] has applied multiplicative and additive renormalizations to define the renormalized chiral condensates ∆ R q .
whered = 0.0232244 and r 1 = 0.3106 fm [6,81]. The chiral condensates at zero temperature ψ ψ q,0 can be determined from ψ ψ q,T at sufficiently low temperatures. An alternative way is to consider them free parameters, determined by fitting lattice data. We find that both yield the same ψ ψ q,0 . Our holographic results are presented in Fig. 9 with a 0 = 30, a 1 = 3, and a l = 0.595 (for light quarks) and a s = 1.23 (for s quark). In particular, note that both condensates change as the temperature increases. Remarkably, the holographic condensates are in quantitative agreement with those given in the lattice QCD simulation [6]. It gives strong support for our holographic setup.
Moreover, let's consider another interesting and nontrivial check by considering the QCD trace anomaly. To do that, we should first convert the chiral condensates ∆ R q to the quark condensates qq T , which can be done by [82] , For QCD with N f = 2 + 1 flavors, the vacuum quark condensates at T = 0 take l l 0 ≡ 0|ll|0 = −[272(5)MeV] 3 for light u, d quarks and ss 0 ≡ 0|ss|0 = −[296 (11)MeV] 3 for s quark [82]. We stress that the quark condensates qq T are different from the chiral condensates ψ ψ q,T . The quark condensates are used to evaluate the QCD trace anomaly. After choosing ∆ R l (∞) ≈ ∆ R l (300MeV) = −0.002 by following lattice QCD [82], we obtain the temperature dependence of quark condensates as shown in Fig. 10. Once again, our holographic results quantitatively agree with the lattice QCD data [6].
The relation between the QCD trace anomaly θ(T ) ≡ (T ) − 3p(T ) and the quark and gluon condensates are given by [83,84] δ β(g) 2g where δf (T ) ≡ f (T ) − f (0) denotes the vacuum subtracted value of the quantity f . The β-function in (38) up to four-loop order reads [85] β(g) 2g = −(2πβ 0 α s + 8π 2 β 1 α 2 s + 32π 3 β 2 α 3 s + 128π 4 β 3 α 4 s ) , where the four-loop QCD running coupling α s in the function of the energy scale µ is given by [86] The corresponding four-loop running quark masses read [87]    We have already known the temperature dependence of the quark condensates qq T using the probe scenario (see Fig. 10) and the trace anomaly θ from the EOS (see Fig. 5). Thus, once knowing the relation between the renormalization scale µ and the temperature T , we can obtain the gluon condensate as a function of T from (38) using our holographic data. Considering that the temperature plays a natural role in a system's energy scale, it is reasonable to identify T as the energy scale µ of (41). More precisely, in our holographic setup, we choose µ/Λ QCD = T /MeV in (40). In Fig. 11, the gluon condensate δ β(g) 2g G 2 T from our holographic setup is denoted by the red dashed curve within the pink band given by the lattice simulation [6]. Moreover, with the four-loop renormalization effect of quark masses (41), the individual contributions from the light and heavy quarks to the QCD trace anomaly quantitatively agree with the lattice data, as shown in Fig. 11. Remarkably, such a simple holographic setup is able to capture the realist QCD thermodynamics and the condensates.
Scenario with first order cosmological QCD phase transition-The CEP predicted by our hQCD model is located at T C = 105MeV, µ C = 555MeV. Therefore, it suggests that the strong first-order phase transition in the cosmic QCD epoch could happen for sufficiently large baryon chemical potential. However, measurements of primordial element abundances and CMB have imposed a strong constraint on the baryon asymmetry η ob B ≈ 10 −10 [69]. The baryon asymmetry is conserved except during baryogenesis and a first-order phase transition. Therefore, to apply to the cosmological QCD . Quark and gluon contributions to the QCD trace anomaly. For the lattice QCD data [6], the pink band denotes the gluon condensate (the right-hand side of (38)), and the solid blue and red curves are, respectively, the contribution from the light (l = u, d) quarks and the s quark to the same equation. The corresponding results from our holographic scenario are shown as red dashed curves.
phase transition, the baryon asymmetry after the firstorder transition in our model should be consistent with η ob B . The value of baryon asymmetry η B at the QCD firstorder phase transition line from our model is presented in Fig. 12. While the value of η B in the false vacuum (green curve) slightly increases as µ B /T increases, the one in the true vacuum (blue curve) decreases quickly. We point out that the constraint on η B by primordial element abundances and CMB only applies to later stages after a first-order phase transition. Thus, we need to compare our η B in the true vacuum at the first order transition to η ob B . Near the CEP, our model gives η B ∼ 10 −2 which is eight orders of magnitude larger than η ob B . Nevertheless, around µ B = 1000MeV (the red square in Fig. 12), η B ≈ 10 −10 reaches the same order of magnitude as the cosmological observation.
Therefore, the cosmological QCD phase transition happens near µ B = 1000MeV (µ B /T ≈ 20) in the first order and agrees with the baryon asymmetry observation. The GWs produced during the QCD phase transition at µ B = 1000MeV are shown in Fig. 4 in the main text. We find that it could be detected via pulsar timing in the future, such as IPTA and SKA.
Moreover, for µ B < 1000MeV, there is a simple scenario with the cosmological QCD phase transition being first order without violating the constraint of a small baryon asymmetry in the later evolution of the universe. In the so-called "little inflation" [70,71] (see also [72][73][74] for earlier work), the universe begins with a large baryon asymmetry and thus crosses the first-order phase transition line. Meanwhile, the large baryon asymmetry is diluted by a short period of inflation to today's observed value. To get an order of magnitude estimate, we approx- imate n Bi (n Bf ) by n B+ (n ob B ) at given µ B with i and f referring to the values before and after the inflation, respectively. We immediately obtain that where a is the scale factor, we have considered that the baryon number in a comoving volume is conserved. In Fig. 13, we show the corresponding e-fold number N ≡ ln(a f /a i ) that is sufficient for the dilution of the baryon number in the little inflation scenario. We point out that this leads to an overestimation of N , since n B has a sudden reduction from the false to true vacuums at the first-order transition. The e-fold number decreases by increasing µ B . Only N ∼ 7 e-folds are sufficient. Nevertheless, the early state of the universe has a relative large baryon asymmetry η B ∼ O(10 −1 ) (see the green curve of Fig. 12). A natural mechanism for generating such a high initial baryon asymmetry is the wellestablished Affleck-Dine baryogenesis [75]. The upper limit for the Affleck-Dine baryogenesis was argued to be η B ∼ O(1) [72]. The case in our model is well below this upper limit. Another scenario that allows a firstorder QCD transition in agreement with observation is considering a large lepton asymmetry [76]. A concrete realization for this scenario was recently given in [77].