Heavy Holographic QCD

We study the phase structure of the Witten-Sakai-Sugimoto model in the plane of temperature and baryon chemical potential, including the effect of a nonzero current quark mass. Our study is performed in the decompactified limit of the model, which, at least regarding the chiral phase transition, appears to be closer to real-world QCD than the original version. Following earlier studies, we account for the quark mass in an effective way based on an open Wilson line operator whose expectation value is identified with the chiral condensate. We find that the quark mass stabilizes a configuration with string sources and point out that this phase plays an important role in the phase diagram. Furthermore, we show that the quark mass breaks up the first-order chiral phase transition curve and introduces critical points to the phase diagram. Similarities of the phase structure to other holographic approaches and to lattice simulations of"heavy QCD"are found and discussed. By making holographic QCD more realistic, our results open the door to a better understanding of real-world strongly coupled hot and dense matter.


Introduction
The Witten-Sakai-Sugimoto model [1,2], which is based on the gauge/gravity duality [3] and in particular on Witten's non-supersymmetric model for low-energy Quantum Chromodynamics (QCD) [4], has been successfully applied to study spectra and properties of mesons, baryons, and glueballs [1,2,[5][6][7]. It has turned out that the model is not only suited for the physics of single particles in the vacuum, but can also be employed to describe thermodynamic systems and study phase transitions. In its original form, the model is dual to QCD at a large number of colors N c (at least in a certain limit, which however is inaccessible with current methods). This is borne out in its phase structure [8], which can be mapped to large-N c QCD predictions [9,10]. Approaching the realistic regime of small N c in a rigorous way is very difficult, although N f /N c corrections, where N f is the number of flavors, have been discussed [11,12]. A promising alternative is the so-called decompactified limit, where the gluon dynamics are decoupled. In contrast to the original form of the model, this limit shows a nontrivial phase structure regarding the chiral phase transition at nonzero chemical potentials, where it behaves not unlike the field-theoretical Nambu-Jona-Lasinio (NJL) model [13]. Examples for applications of the Witten-Sakai-Sugimoto model in that limit to phases of QCD are the discovery of "inverse magnetic catalysis" in dense matter [14], and the possibility of a quark-hadron continuity at large chemical potential [15]. These studies have all been performed in the chiral limit, i.e., at vanishing current quark mass (a constituent quark mass is induced through spontaneous chiral symmetry breaking). In this paper, our goal is to go beyond the chiral limit and study the effect of an explicit chiral symmetry breaking on the phase diagram of the model in the decompactified limit.
The Witten-Sakai-Sugimoto model is based on a setup of N c D4-branes and N f pairs of D8-and D8-branes. These "flavor branes" represent left-and right-handed fermions and, for N c N f , can be treated as "probe branes", i.e., they do not backreact on the background geometry determined by the D4-branes. Including a current quark mass in this setup is not straightforward because there are no transverse directions that provide a massgap for excitations of strings that have endpoints at one color and one flavor brane. Two approaches have been proposed to circumvent this problem, both relying on the physics of open strings connecting the D8-and D8-branes. In the first approach a bi-fundamental scalar field is introduced, which becomes tachyonic in the infrared and which is used to account for the chiral condensate and the current quark mass [16][17][18]. We follow the second approach, which is based on an open Wilson line between the left-and right-handed flavor branes [19][20][21][22][23][24]. The expectation value of this operator is comparable -but due to its nonlocal nature not identical -to the usual chiral condensate. While this expectation value was calculated and discussed in previous works, it has, to the best of our knowledge, not yet been implemented into a fully consistent calculation of the phase structure, which includes its backreaction on the embedding of the flavor branes. Here we provide such a calculation for nonzero temperatures and chemical potentials. We emphasize that this treatment of the current quark mass is not a rigorous top-down approach. In analogy to ordinary chiral perturbation theory we simply add a term to the action that has the form quark mass times chiral condensate. We shall discuss the novelties in the phase structure of the system together with the conceptual problems that arise in this approach, keeping in mind that results for very large current quark masses have to be treated with care since higher order mass terms are neglected from the beginning.
Our work is related to various previous studies. Most notably in the holographic context, the phase structure of strongly interacting matter has been studied in the D3-D7 setup [25], where a current quark mass is straightforwardly included and the interpretation of the chiral condensate is unambiguous [26][27][28]. Despite the differences in how to introduce the quark mass we shall see that both from the holographic point of view as well as in the resulting phase diagrams we find intriguing parallels with the D3-D7 system, which are less striking in the massless limit. Another holographic approach used to address similar questions is a bottom-up approach that works in the Veneziano limit, where the ratio N f /N c is held fixed while sending the number of colors to infinity [29,30]. For example, we shall compute the speed of sound for nonzero temperatures and chemical potentials, which also has been done in that approach [31]. Last but not least, we shall compare our results to recent calculations of QCD on the lattice, where a combined expansion for large coupling and large quark masses is used to circumvent the sign problem and thus is able to provide results for large baryon chemical potentials [10,32].
Our paper is organized as follows. In Sec. 2 we introduce the model and discuss the mass correction to the action. The solutions to the equations of motion and to the stationarity equations of the free energy density are derived in Sec. 3 for three different embeddings of the flavor branes. The numerical evaluation with all physical results is presented and discussed in Sec. 4. Or main result is given in Sec. 4.3 in the form of phase diagrams for different values of the current quark mass. In Sec. 5 we summarize and give an outlook.

Brief introduction to the model and Dirac-Born-Infeld action
We start with a very brief description of the model, to establish the notation and the action to which we add the quark mass correction. More detailed accounts of the Witten-Sakai-Sugimoto model can be found in the original works [1,2] and in the review [33]. The model is based on a type-IIA string theory background generated by a large number N c of D4-branes extended over the directions (x 0 , x, x 4 ), where (x 0 , x) are the 4d space-time directions. As proposed by Witten [4], the x 4 direction is compactified on a Kaluza-Klein (KK) circle of radius M −1 KK with boundary conditions for the fermions that lead to the breaking of supersymmetry. As we are interested in thermal properties of the system, also the imaginary time direction τ = ix 0 is compactified, with radius (2πT ) −1 , where T is the temperature. Thus, while the x 4 radius is a fixed dimensionful parameter of the model, the τ radius becomes smaller with increasing temperature. The AdS/CFT duality connects the low-energy description of this system with a non-supersymmetric and non-conformal Yang-Mills field theory. The parameters of the string theory are related to the field-theoretical ones via where R is the curvature radius of the bulk geometry, s is the string length, g s is the string coupling, and λ is the 't Hooft coupling. A supergravity description is accurate for small curvature, i.e., for large values of λ. At low temperatures, the topology of the manifold spanned by x 4 and the radial ("holographic") coordinate u is that of a cigar, while the manifold spanned by τ and u is cylindershaped. Above the critical temperature T c = M KK /(2π) the energetically favored background geometry changes and τ and x 4 exchange their roles: a black hole appears. This Hawking-Page-type transition on the gravity side is usually identified with the deconfinement transition of the gauge theory that lives at the holographic boundary u = ∞ (see however [34]), and we shall refer to the background where the τ -u-manifold is cigar-shaped as the "deconfined geometry". (Strictly speaking, for temperatures larger than the Kaluza-Klein scale, unwanted modes become dynamical and the effective four-dimensional description breaks down.) In the deconfined geometry, the topology of the manifold spanned by x 4 and u is a cylinder, as illustrated in Fig. 1. This figure also includes the flavor branes (to be discussed momentarily) and their three different stable embeddings in the presence of a quark mass (to be discussed in Sec. 3). We will always work in the deconfined background, The asymptotic distance of the flavor branes is given by the dimensionless parameter (such that the x 4 radius is 1 and = π corresponds to an antipodal embedding). In the mesonic and LTQ phases the flavor branes connect above the horizon u T at a point that we denote by u c and which has to be determined dynamically.
The mesonic phase has a smooth embedding and vanishing baryon number. In the LTQ phase, which has nonzero baryon number, there is a cusp at which strings are attached, here represented by a straight line. In the HTQ phase, which also has nonzero baryon number, the flavor branes extend all the way to the horizon. The dashed lines indicate the massless limit m q = 0 in which there is no stable string configuration and in which the stable phase with nonzero baryon number has straight, disconnected flavor branes. whose metric is where dΩ 2 4 is the metric of a unit 4-sphere, and the blackening function f T is where u T is the position of the horizon in the black hole geometry, related to the dimensionless temperature t via Here and throughout the paper we work with the same dimensionless quantities used before in a series of works [15,35,36]. The relation to their dimensionful counterparts involves powers of M KK and the curvature radius R as well as some numerical constants, see for instance Table I in Ref. [35]. Here we restrict ourselves to giving the definitions only of the most important quantities, which are needed to translate our results into physical units. This is done in Table 1, where we have introduced the useful abbreviation temperature quark chemical quark number constituent f 2 π m 2 π potential density quark mass Table 1. How to obtain the physical, dimensionful quantities from their dimensionless counterparts t, µ, n, M q , α that are used in the calculation. In all plots in Sec. 4 we use an additional rescaling with the asymptotic separation of the flavor branes,t = t ,μ = µ 2 ,ñ = n 5 ,Ω = Ω 7 ,M q = M q 2 , α = α 4 ,λ = λ/ .
As usual for holographic models, fundamental matter is introduced by adding probe D-branes. Sakai and Sugimoto included N f D8-D8 pairs located at the antipodes of the x 4 circle. The U (N f ) L × U (N f ) R gauge symmetry on the flavor branes gives rise to chiral symmetry in the field theory. A quark chemical potential is included as the boundary value of the temporal component of the abelian gauge field, such that the baryon chemical potential is N c µ (we include a hat in the notation of the abelian gauge field for notational consistency with previous works). Below the critical temperature for deconfinement, the topology of the background forces the D8-D8 pair to join in the bulk, forming a ∪-shaped configuration. This provides a geometric realization of the spontaneous breaking of chiral symmetry to the diagonal U (N f ). We shall refer to this phase as "mesonic" since the baryon density is zero. Baryons can be added in the form of instantons in the ∪-shaped configuration, but here we shall not include this possibility. Above the critical temperature, the probe branes are straight and end at the horizon. Since they are disconnected, gauge transformations on D8-and D8-branes are independent and chiral symmetry remains intact. This phase does have a nonzero baryon density, generated by deconfined quarks. Chiral restoration and deconfinement occur at the same µ-independent critical temperature in the case of an antipodal separation of the flavor branes. However, this separation L can be considered as a second geometric parameter of the model, with L = π/M KK corresponding to antipodal separation. We shall work with the dimensionless version = LM KK , i.e., where x 4 (u) describes the shape of (half of) the embedding of the probe branes. Therefore, the parameters of our model are λ, M KK and . Our main results are of qualitative nature, and thus we do not have to choose specific values for these parameters. Only towards the end of Sec. 4.3 we shall briefly discuss some quantitative results using a specific parameter set. Non-antipodal asymptotic separations, < π, generalize the model [37] and lead to interesting new results. In particular, for sufficiently small , namely < 0.30768π, the physics of chiral symmetry breaking and confinement get disentangled because the ∪-shaped embedding now becomes a solution also in the deconfined geometry, and a free energy comparison between the different embeddings (straight and connected) shows that the chiral phase transition depends on the chemical potential, as expected in real-world N c = 3 QCD. We work in the limit π, where the analysis of chiral symmetry breaking is effectively decoupled from confinement effects, and the model is not unlike a holographic version of a NJL model [13,38,39]. In particular, we shall work with the deconfined geometry (2.2) for all temperatures, having in mind that by decreasing L at given M KK or increasing the x 4 radius M −1 KK at given L ("decompactification") we can tune the confinement scale to arbitrarily small values compared to the scale for chiral symmetry breaking. Of course, by taking this limit we lose some of the top-down control of the original model since now Kaluza-Klein modes become potentially relevant. The gain, however, is a much richer phase structure which has proven to yield interesting physics that are, at the very least, comparable to well-established field-theoretical models. It has been shown, for instance, that the parameters of the model in the decompactified limit can be fitted to reproduce properties of nuclear matter at saturation density [15], making it a viable candidate for the study of dense matter in neutron stars.
In the chiral limit, the different probe brane embeddings and their free energies are obtained from the (Euclidean) Dirac-Born-Infeld (DBI) action which is derived from the induced metric on the D8-branes and the dilaton field. Here, primes denote derivatives with respect to u, and we have abbreviated where α = 2 s , T D8 = (2π) −8 −9 s is the D8-brane tension and V S 4 = 8π 2 /3 is the volume of the unit four-sphere, while the factor 2 in the first expression accounts for the two halves of the brane worldvolume. If the flavor branes reach the horizon, the lower boundary of the integral in Eq. (2.8) has to be replaced by u T . We have used translation symmetry in the flat spacetime directions to extract the overall factors of the 3-volume V and inverse temperature. Due to the flavor symmetry, the number of flavors simply appears as a prefactor in front of the action. The grand-canonical potential is given by which defines the dimensionless free energy density Ω as a function of the dimensionless thermodynamic variables µ and t. The resulting phase diagram in the µ-t-plane shows a first-order chiral phase transition [40], separating the chirally broken, ∪-shaped mesonic phase from the chirally restored phase with straight and disconnected branes. We shall discuss how this result is altered by a nonzero current quark mass, which we introduce now.

Quark mass correction
We first recall that, while fundamental quarks are represented by strings stretched between the D4-and the D8-branes, mesons can be described by strings with both endpoints on the D8-branes, and they manifest themselves as scalar and vector fluctuations in the worldvolume theory. When the D8-D8 pairs join in the interior, one finds massless scalar modes in the spectrum. These modes are the Goldstone modes generated by the spontaneous breaking of chiral symmetry. Their fluctuations can be studied in terms of the chiral field , where x stands for the spacetime coordinates and f π is the pion decay constant. For instance, for N f = 2, π(x) accounts for the three massless pion modes (and the η, which however acquires a mass through the axial anomaly). In the holographic theory the chiral field corresponds to the holonomy where u = (u 3 c + u c z 2 ) 1/3 , such that the coordinate z parameterizes the entire connected flavor branes, where a z is the (dimensionless) radial component of the non-abelian gauge field, and where P denotes path ordering.
In QCD, the pions are of course not massless because chiral symmetry is explicitly broken. In the Witten-Sakai-Sugimoto model, explicit chiral symmetry breaking is not as straightforward as in other holographic approaches such as the D3-D7 setup. In the latter, one only needs to impose a fixed asymptotic separation between flavor and color branes along one of the transversal directions. The radial dependence of this D3-D7 distance is captured by one of the scalar fields of the worldvolume theory. As usual, once the equations of motion are solved, the corresponding wavefunction allows one to identify a nonnormalizable and a normalizable mode, whose coefficients are proportional to the current quark mass m q and to the expectation value of the quark bilinear, i.e., the chiral condensate qq , respectively. In the Witten-Sakai-Sugimoto model, however, there is no transverse direction to separate the color from the flavor branes. Moreover, left-and right-handed fermions are separated by the extra dimension x 4 , such that there is no obvious local expression for the chiral condensate. Following Ref. [19], we consider instead the non-local, but gauge-invariant, operator given by the open Wilson line where i, j = 1, . . . , N f are flavor indices, and consider the following correction to the action, where the trace is taken over flavor space, and we have assumed all quark flavors to have the same mass m q . This term has the same form as the lowest-order mass term in chiral perturbation theory where M is the mass matrix, M = diag(m u , m d ) for N f = 2, and U is the chiral field, represented by Eq. (2.11) in the holographic theory. The prefactor qq 0 is the chiral condensate in the vacuum at vanishing quark mass. With m q = (m u + m d )/2 it is related to the pion decay constant f π and the pion mass m π in the vacuum via the Gell-Mann-Oakes-Renner relation [41] − m q qq 0 = f 2 π m 2 π 2 . (2.15) The mass term (2.14) is also used in the approach mentioned in the introduction based on the tachyonic field and moreover was employed (independent of the approach for the chiral condensate) in the study of glueball decay rates in the Witten-Sakai-Sugimoto model [42]. The connection between (2.13) and (2.14) is made as follows. At strong coupling, the holographic dual of Eq. (2.13) can be written in terms of the semi-classical expression for the expectation value O . As usual, the open Wilson line corresponds to the Euclidean worldsheet of a fundamental string which ends at the boundary line, and its one-point function is proportional to the exponential of the area of the associated minimal surface, O ∝ e −Sws . Including boundary terms, we write the worldsheet action as S ws = S NG +S ∂ , where S NG is the Nambu-Goto action, and the boundary contribution S ∂ = i ∞ −∞ dz a z gives rise to the chiral field (2.11) and thus encodes the meson fluctuations. A nontrivial a z also arises in the presence of baryons, which are instantons on the flavor branes. Therefore, this setup can also be used to compute quark mass corrections to the baryon spectrum [23,43,44]. We identify the chiral condensate with the Nambu-Goto part, O = − qq U , where qq = −ce −S NG with a proportionality constant c (of mass dimension 3). The minus sign is added for convenience such that c is positive. The chiral condensate in general contains medium and mass corrections and is related to the vacuum chiral condensate at where S 0 NG is the Nambu-Goto action evaluated in the vacuum and at zero quark mass, which has a simple analytical form, see Eq. (3.17). We shall use Eq. (2.16) later to compute the normalized chiral condensate. Using Eqs. (2.9), (2.15), and assuming the chiral condensate to be uniform, we can write Eq. (2.13) as where, in the second step, we have neglected the mesonic fluctuations, i.e., U = 1, and performed the (then trivial) trace over flavor space and the spacetime integral. We have introduced the dimensionless prefactor α, which is only nonzero in the presence of a current quark mass m q , In other words, if we write the chiral field U in Eq. (2.17) in terms of the pion field π(x), then the prefactor of the quadratic term in π(x) gives the square of the pion mass. The resulting expression for m 2 π is nothing but Eq. (2.18). Despite the explicit appearance of N c , the "mass parameter" α does not scale with N c because O and thus c are proportional to N c [21]. Notice also that we do not use any explicit form of c. As a consequence, our mass parameter cannot be mapped directly onto the current quark mass m q , but only to the combination f 2 π m 2 π ∝ N c . We will study how the mass perturbation (2.17) affects the embedding of the flavor branes for different values of the mass parameter (2.18), while neglecting its effect on the background geometry, which is suppressed by a power of g s [19].
It remains to compute the Nambu-Goto action. We have where the factor 2 accounts for the two halves of the connected flavor branes, and the prefactor λ 0 arises from combining the prefactor (2πα ) −1 that appears in the Nambu-Goto action with the factors from our use of dimensionless coordinates. With the metric (2.2) and introducing an ultraviolet cutoff Λ, Eq. (2.19) becomes where we have integrated by parts, used the boundary condition (2.7), and introduced where 2 F 1 denotes the hypergeometric function. For the boundary term we have used φ T (u) u for u u T . At zero temperature, φ T (u) = u for all u. Later we will also need In Eq. (2.20) we have allowed for a nonzero x 4 (u c ), which allows for the flavor branes to be connected with a straight segment constant in u. We shall see that this becomes relevant for the case where the flavor branes approach the horizon, see Sec. 3.3. The action (2.20) diverges for Λ → ∞. We can simply renormalize it by subtracting the vacuum contribution, i.e., by dropping the contribution proportional to Λ, which does not depend on temperature or chemical potential. In Ref. [24] a temperature-dependent renormalization was employed, which amounts to subtracting the Nambu-Goto action for disconnected, straight branes, where x 4 = 0 and u c has to be replaced by u T , Since a theory should be renormalized in the vacuum, before switching on any medium effects, this seems incorrect from a theoretical point of view, as argued for instance in Ref. [45] in a similar context, see also footnote 6 in Ref. [46]. The reason why one might be tempted to subtract S || NG is that in this case the renormalized worldsheet action is proportional to the area between the horizon and the holographic boundary that is complementary to the area enclosed by the flavor branes. It is then conceivable to interpret this area as an order parameter for chiral symmetry breaking. However, even in this approach the order parameter is not proportional to the chiral condensate because the worldsheet action appears in the exponent. (For that reason it was even suggested to interpret O − 1 as an order parameter [24], and also the possibility of using ln O was mentioned in a footnote in Ref. [19].) Here we follow the theoretically more straightforward way to subtract the vacuum contribution in Eq. (2.19), keeping in mind that the interpretation of the chiral condensate has to be taken with care. This might not be too surprising since it is constructed from an open Wilson line operator, which is not identical to a local quark-antiquark operator to begin with. We shall come back to this point when we present the numerical results, see Sec. 4.1 and in particular Fig. 3.
To summarize, and to write the mass term in a convenient form for the following calculation, we have arrived at with (after renormalization) where the u-integral is finite and we have thus reinstated Λ → ∞. While A is constant in the holographic coordinate u, it is a functional of x 4 (u) and thus implicitly depends on temperature and chemical potential (it also depends explicitly on temperature). We shall see that via the equations of motion x 4 (u) also depends on A itself, such that Eq.
(2.25) is an implicit equation for A, which we have to solve numerically. The mass term (2.24) introduces -in addition to the obvious dependence on the mass parameter α -also a nontrivial dependence on the 't Hooft coupling. In the absence of a quark mass, λ can be formally eliminated by an appropriate rescaling of the variables. This is no longer true in the massive case, and we thus have to discuss the phase structure of the model in the two-dimensional parameter space spanned by λ and α.

Different probe brane configurations
In this section we discuss three different probe brane configurations, which appear as classical solutions of the action Here, in addition to the DBI action and the mass correction discussed above, we have added a source contribution, first considered in Ref. [8]. The sources are given by static fundamental strings stretching from the tip of the flavor branes all the way down to the horizon (and distributed uniformly along the 3d flat directions), see Fig. 1 (middle). The endpoints of these strings at u = u c act as charged sources for the gauge potentialâ 0 (u), therefore their number density can be identified with the quark number density 2 . In our conventions, the number of strings N s is related to the dimensionless baryon number density n by and we have Here, the first term is the Nambu-Goto contribution (u c −u T simply being the string length), while the second term is the boundary contribution 3 . From Eq. (3.3) we can read off the energy of a single string, which corresponds to the constituent quark mass. We may thus define the dimensionless version of the constituent quark mass, in accordance with previous works [19,22], as Collecting the constants in Eq. (3.3) then gives the dimensionful version, see Table 1.
The equations of motion for the abelian gauge field and the embedding function, become in integrated form with A from Eq. (2.25) and k being an integration constant. In the presence of string sources, n appears directly due to (3.3); when they are absent, n arises as an integration constant. We can solve Eqs. (3.6) algebraically for x 4 andâ 0 , where we have abbreviated This yields the asymptotic expansions at the holographic boundarŷ (3.9b) Via the AdS/CFT dictionary we see that n is indeed the number density associated with the chemical potential µ. Note that the quark mass changes the asymptotic behavior of the embedding function x 4 (u) since in the presence of a quark mass A > 0. The boundary conditions at u = ∞ are given by Eqs. (2.6) and (2.7). For the following it is useful to write them in the form Finally, the dimensionless free energy density, defined in Eq. (2.10), is obtained by inserting the solutions (3.7) back into the Euclidean action. This yields The free energy depends on the thermodynamic variables µ and t, on the externally given parameters , λ 0 , α (the dependence on M KK is absorbed in the definition of the dimensionless quantities), and on the variables k, n, u c ,â 0 (u c ), x 4 (u c ), which have to be determined dynamically, depending on the particular boundary conditions we choose at u = u c . We shall now discuss separately the three phases schematically shown in Fig. 1.

∪-shaped embedding: mesonic phase
First, we consider a ∪-shaped embedding, where the branes connect smoothly at u = u c , i.e., we have x 4 (u c ) = 0 and x 4 (u c ) = ∞. Then, assuming a finiteâ 0 (u c ), Eq. (3.6a) implies n = 0, which, in turn, with Eq. (3.7a), gives a constant abelian gauge field,â 0 (u) = µ. For this solution, the string sources play no role, and so far there is no qualitative difference to the massless case. We can rewrite the free energy (3.11) as where we have used Eqs. (3.7b), (3.10b), and n = 0. This free energy does not depend on µ, which is obvious on physical grounds since the baryon density is zero. The only way to introduce baryon number in this phase is through actual baryons, which we ignore throughout this paper. The form (3.12) is particularly useful for finding the stationary points of Ω with respect to the remaining variables k and u c . We find that the stationarity condition with respect to k is equivalent to the boundary condition (3.10b), while the stationarity condition with respect to u c reads . (3.13) Since the square root appears also in the denominator of x 4 , see Eq. (3.7b), stationarity with respect to u c turns out to be equivalent to the smoothness condition for the embedding at u = u c , i.e., x 4 (u c ) = ∞. We use this condition to express k in terms of u c , This reduces the evaluation of the mesonic phase to numerically solving for A and u c and inserting the solution into the free energy (3.12).
Here and throughout the paper we do not expand our equations for small α. Such an expansion would be consistent with our starting point, since our quark mass correction to the action is linear in m q . However, we have checked that not much is gained from this expansion, the resulting equations are even somewhat more tedious and do not facilitate the numerical evaluation. Therefore, we formally keep all orders in the mass, but need to keep in mind that results for very large quark masses can only be obtained as an extrapolation of our approximation.
An analytical result can be obtained from Eqs. (3.15) in the massless limit at zero temperature. With A = 0, f T (u) = 1, and φ T (u) = u we have with the new integration variable u = u/u c (and dropping the prime again) where, in the second equation, we have used the result for u c from the first equation. Consequently, in the vacuum and at zero quark mass we have This result is identical to that of the confined geometry with non-antipodal flavor branes, see Eq. (4.20) of Ref. [21], where this result is expressed in terms of the five-dimensional 't Hooft coupling λ 5 = 2πλ/M KK , and where subleading corrections from the Fradkin-Tseytlin contribution were computed additionally (including the factor N c mentioned below Eq. (2.18)).

3.2
-shaped embedding: low-temperature quark phase Next we consider a -shaped embedding in which the flavor branes connect in a cusp that is caused by the string sources. This phase has nonzero baryon number, but since there are no baryons in the system, baryon number must be created from quarks. Moreover, we shall see that this phase does not exist at large temperatures. Hence we refer to it as the low-temperature quark (LTQ) phase, see Fig. 1 (middle). We use the general form of x 4 andâ 0 from Eqs. (3.7), as well as Eqs. (3.10) and x 4 (u c ) = 0 to write the free energy (3.11) as Again, we need to consider the stationarity equations of the free energy, this time with respect to k, u c , and n. Stationarity with respect to all three variables ensures thermodynamic consistency: the derivative of Ω = Ω [k(µ, t), u c (µ, t), n(µ, t), µ, t] with respect to µ at fixed t is identical to the explicit derivative if all derivatives with respect to k, u c , n are zero. And, as Eq. (3.18) shows, the explicit derivative is −n. Therefore, the identification of n as the baryon density from the AdS/CFT dictionary is in accordance with its thermodynamic definition, as it should be. We find that stationarity with respect to k is equivalent to Eq. (3.10b), as in the mesonic phase, and thus does not yield any new information. The stationarity equations with respect to n and u c read (3.19b) The first equation obviously impliesâ 20) which means that S source does not yield a contribution to the free energy (although it affects the solution indirectly through the equations of motion). From the second relation (3.19b) we derive a condition for k which takes exactly the same form as in the mesonic case, see Eq. (3.14). In particular, the baryon density n does not enter this relation. The only condition for n from Eq. (3.19b) is n > 0, meaning that we always deal with a net baryon, not anti-baryon, density. Inserting the result for k into Eq. (3.7b) yields which confirms that there is a cusp at the tip of the connected flavor branes, x 4 (u c ) < ∞, created by the string sources.
Having eliminated k, we are left with solving for A, u c , and n for given µ and t and inserting the results back into the free energy (3.18).
In the limit n → 0, the first two relations are identical to those of the mesonic phase and, as Eq. (3.21) shows, the brane embedding becomes smooth. This limit is assumed at a certain value of µ which is given by the third relation: since at zero density we haveâ 0 (u) = 0, the LTQ phase approaches the mesonic phase at the point µ = u c − u T , where u c of course depends implicitly on µ. In other words, at this particular µ it becomes possible to populate the system with quarks. This confirms the interpretation of u c −u T as the constituent quark mass, as defined in Eq. (3.4). The continuous geometrical connection between the mesonic and LTQ phases suggests that a continuous phase transition between them is possible. We shall see later that such a transition is realized for sufficiently large quark masses. The string solution discussed here was already considered in the massless case [8]. It was argued that this phase is unstable due to a negative number susceptibility. Our results agree with those of this reference, but we have found that even in the massless limit there is a regime in which the LTQ phase has positive number susceptibility: if the number density in Fig. 11 of Ref. [8] is continued to lower values of µ one finds a stable branch. Nevertheless, in the absence of a current quark mass, the LTQ phase is metastable at best, i.e., even when it has positive number susceptibility it is energetically disfavored. This is different in the presence of a quark mass, and we shall see that the LTQ phase plays a key role in the phase structure of the model. In fact, we will see that this phase behaves very similarly to a certain solution in the D3-D7 approach [25,26,47], where the flavor branes connect to the horizon in a "long spike" that "resembles a bundle of strings" [25]. The analogy exists not only in the bulk geometry, but will also become manifest by comparing our phase diagrams with the corresponding ones in the D3-D7 model, see Sec. 4.3. This observation suggests that the string solution considered here might be an approximation to a more complicated, but smooth, embedding of the D8-D8 pair, where the cusp and the string sources are replaced by a spike-like shape of the flavor branes, see also Refs. [48,49]. We leave it to future studies to identify such a solution in the Witten-Sakai-Sugimoto model.

3.3
-shaped embedding: high-temperature quark phase Finally, we consider a configuration that can be understood as a deformation of the straightbrane configuration in the chiral limit, see Fig. 1 (right). This phase also has nonzero baryon number, which is created by quarks. We shall see that this solution exists only for sufficiently large temperatures, and thus we refer to it as the high-temperature quark (HTQ) phase 4 . For m q = 0 the brane embedding at non-zero quark density is straight, i.e., in this case x 4 (u) = 0 due to k = A = 0, see Eq. (3.7b). Now, for m q > 0 and thus A > 0, x 4 (u) = 0 is not a solution anymore. We still find a solution, however, where the branes extend all the way down to the horizon. Therefore, for the HTQ phase, we set u c = u T and consider x 4 (u T ) as a variable that has to be determined dynamically (for straight branes, x 4 (u T ) = /2). Then, the (renormalized) Nambu-Goto action reads where now, as opposed to the previous cases, the lower boundary term is nonzero, with φ T (u T ) being proportional to u T , see Eq. (2.22). As for the straight-brane solution we requireâ 0 (u T ) = 0. Rewriting the free energy (3.11) with the help of Eqs. (3.10) and replacing the lower boundary of the radial integration by u T yields Now the stationarity conditions with respect to n and k are equivalent to the boundary conditions (3.10a) and (3.10b), respectively. Hence the only relevant stationarity equation where the second term originates from the explicit dependence of A on x 4 (u T ). We conclude that k = −Aφ T (u T ) .

(3.26)
Near the horizon we find i.e., the curved flavor branes are horizontal at the horizon, as indicated in Fig. 1. This is in contrast to the massless case, where the straight branes approach the horizon vertically.
We are left with solving (3.28) for A and n, while x 4 (u T ) can be computed from (3.10b) afterwards if needed. Of course, only solutions with x 4 (u T ) > 0 make sense, otherwise the branes would cross each other. We find that indeed for small temperatures one runs into solutions with x 4 (u T ) < 0, which have to be discarded. As a result, there is a (µ-dependent) temperature below which the HTQ phase does not exist. This is the temperature at which the -shaped HTQ configuration approaches the -shaped LTQ solution. At this point x 4 (u T ) = 0, and it is easy to see that all the conditions we have derived are nothing but the u c → u T limit of those described in the previous subsection. We have thus shown that, on the one hand, the LTQ solution connects continuously to the mesonic solution for n → 0 and, on the other hand, the HTQ solution connects continuously to the LTQ solution for x 4 (u T ) → 0. In other words, with the -shaped string solution we can interpolate continuously between the ∪-and -shaped solutions. This geometric continuity suggests that continuous phase transitions mesonic-LTQ and LTQ-HTQ are conceivable. Whether these transitions are smooth in terms of higher derivatives of the free energy has to be checked numerically, and we will discuss this point in the next section.
One might ask in which sense the HTQ solution breaks chiral symmetry explicitly. After all, we have introduced a current quark mass, and thus neither the action nor any of the stable states should be chirally symmetric. Geometrically, chiral symmetry in the Witten-Sakai-Sugimoto model is intact (broken) if the flavor branes are disconnected (connected). Now, just as the straight-brane solution in the chiral limit, one might think that the HTQ configuration consists of disconnected branes. However, it might be more accurate to think of the branes in the HTQ phase to be connected by a straight segment at fixed u = u T from x 4 = −x 4 (u T ) to x 4 = +x 4 (u T ) along whichâ 0 (u T ) = 0. This segment gives no contribution to the on-shell action, as can be seen for instance by changing variables to u(x 4 ), and thus the above results are independent of this consideration. The fact that x 4 (u T ) = ∞ further suggests that the branes "want" to connect, at least the connecting piece can be inserted smoothly. Connectedness (although not smoothness) is therefore reflected in our choice of the subscript for the HTQ phase. Another argument for the branes to be connected is the fact that for disconnected branes the expectation value of the open Wilson line operator has been argued to vanish, O = 0 [21,22]. Interpreted as an exactly vanishing chiral condensate, this would suggest chiral symmetry to be exactly restored, which, in the presence of a current quark mass, is only expected for asymptotically large temperatures or chemical potentials.

Results
The remainder of the paper is devoted to the numerical evaluation of the phases just discussed and to the resulting phase diagrams obtained for different values of the quark mass parameter α and the 't Hooft coupling λ. We have reduced the numerical evaluation to solving systems of algebraic equations, which contain numerical integrals. One can easily show that the asymptotic separation of the flavor branes can be eliminated from the equations in each phase by rescaling all quantities with appropriate powers of . As a consequence, we will plot only rescaled quantities, which we denote by a tilde. To facilitate translation to the original, un-rescaled quantities, we have collected the definitions of all relevant rescaled quantities in the caption of Table 1. This table therefore shows how to translate the dimensionless, rescaled variables of all following results into concrete physical values for given model parameters λ, , M KK .
We start by discussing the thermodynamics and phase transitions by fixing either µ or t in Secs. 4.1 and 4.2. Then, we present the phase diagrams in the t-µ-plane in Sec. 4.3. Most results are shown for the value of the rescaled 't Hooft couplingλ ≡ λ/ = 15, focusing on changes that occur upon variation of the rescaled mass parameterα ≡ α 4 . The particular value forλ is chosen somewhat arbitrarily, but it is almost identical to the value chosen in Ref. [15], obtained from fitting the model parameters to properties of nuclear matter at saturation. (This fit will receive corrections, however, if quark mass effects are included.) At the end we shall also discuss how variations inλ affect the phase structure.  In Fig. 2 we see that in the massless case the mesonic solution is two-valued for all temperatures belowt 0.17, above which no solution exists. In this case, the straightbrane solution exists for all temperatures and a first-order transition between the two phases occurs. Here, we have shown this solution forμ = 0. At larger chemical potentials the straight-brane solution, while still existing for all temperatures, moves to lower values of Ω, such that the phase transition moves to lower temperatures until the chirally symmetric phase is preferred for allt ifμ is sufficiently large. This results in the phase diagram of the model in the decompactified limit at zero current quark masses, which is well known [40].

Thermodynamics and phase transitions
For nonzero quark mass, the mesonic solution becomes single-valued for small temperatures since the unstable branch does not reach back all the way to zero temperature, i.e., the point where the tip of the connected flavor branes touches the horizon is now at nonzero temperature. At vanishing chemical potential (left panel), not much is changed compared to the massless scenario: the HTQ phase connects to the unstable branch of the mesonic phase and again there is a first-order transition between the two phases. We have checked that the discontinuity of this transition becomes smaller as we increase the quark mass parameter, but a continuous transition is never obtained. In fact, ifα is chosen too large, the mesonic solution ceases to exist even for small temperatures. This would imply that there is no stable phase in a certain region of the phase diagram, which is unphysical. As discussed above, we cannot trust our approximation to arbitrarily large masses, and thus this behavior is not surprising. We shall therefore only consider mass parametersα -more precisely, only parameter pairs (λ,α) -for which we can determine a stable equilibrium state for allt andμ.
A more interesting situation, qualitatively different from the massless case, occurs for nonzero chemical potentials, as shown in the right panel of Fig. 2. Here, the mesonic and HTQ solutions do not connect continuously anymore. For sufficiently largeμ we find that at the endpoint of the HTQ solution this solution has smaller free energy than the mesonic solution. Without a third solution, this would imply an unphysical jump in the free energy. This third solution is the string solution: it connects with the mesonic one where µ = u c − u T and n = 0, and with the HTQ solution where u c = u T at some nonzero n. As a consequence, the right panel of the figure shows a first order phase transition from the mesonic to the LTQ phase followed by a transition to the HTQ phase. From that figure it is difficult to say how smooth the LTQ-HTQ transition is. Our results strongly suggest that the density, i.e., the first derivative of the free energy, is continuous (see also Fig.  4). Therefore, it is not a first-order transition. We shall later come back to the question whether it is a higher-order transition or a smooth crossover.
In Fig. 3 we plot the chiral condensate (2.16) (using Eq. (3.17)) and the constituent quark mass (3.4) as a function ofμ for two different temperatures and four different values of the mass parameter. The values of the mass parameter are the same as used later for the phase diagrams in Sec. 4.3. We see that at zero temperature there is a transition from the mesonic to the LTQ phase. This transition is of first order for small quark masses and becomes continuous for larger quark masses. The HTQ phase does not exist at zero temperature. At nonzero temperatures, heret = 0.1, the HTQ phase starts playing a role, and we find a mesonic-HTQ first order transition for small quark masses and the sequence mesonic-LTQ-HTQ, as already seen in Fig. 2 as a function oft. The curves shown here do not exhaust all possible scenarios. In particular, for larger temperatures there is a first-order LTQ-HTQ transition. This will become clear in Fig. 4 and in the phase diagrams of the next subsection. In Fig. 3 we also see that both chiral condensate and constituent quark mass (and in fact all thermodynamic quantities) approach the same value at large chemical potential, independent of the value of the quark mass. This is expected of course since for sufficiently large chemical potentials (and also for large temperatures) the quark mass should be negligible. Geometrically speaking, this means that the flavor branes in the HTQ phase become more and more straight and approach the chirally symmetric embedding. In the special caset = 0 it is the LTQ phase that approaches the straight-brane embedding: as µ becomes larger the connected branes approach the horizon while assuming a rectangular shape. These asymptotic limits are a good check for our numerical calculation. All these observations could also have been made by plotting different thermodynamic quantities, for instance the density.
Let us now comment on the quantities themselves that are plotted in Fig. 3. Chiral condensate and constituent quark mass behave qualitatively similar, but are not connected by a simple relation. In the simplest version of the NJL model, they are proportional to each other, connected by a (dimensionful) coupling constant. In QCD, the constituent quark mass does not only originate from the chiral condensate since gluons play an important role as well. It is thus not unphysical to observe different behaviors of the two quantities. However, we have already argued that our "chiral condensate" is really constructed from a non-local operator and its interpretation is thus not obvious. Indeed, we observe for instance that asymptotically, where one might expect chiral symmetry to be effectively restored, we obtain a nonzero value. This value is indicated by a horizontal dashed line and is given by where the arrow stands for the limit µ → ∞ or t → ∞, and where we have used Eqs. (2.4), (2.22), and (3.28). At fixed λ 0 and , this expression grows exponentially with t. This t-dependence is left since we have not renormalized the Nambu-Goto part of the worldsheet action by the result of the straight-brane solution but with a vacuum term that does not depend on temperature, see discussion around Eq. (2.23). Interestingly, however, if is sent to zero at fixed λ 0 and t, the asymptotic chiral condensate, compared to the vacuum value, approaches zero, as one would expect. (Eq. (4.1) is only valid for the LTQ and HTQ phases, i.e., the mesonic piece of the curves in Fig. 3 does not at the same time go to zero.) The difficulties in the interpretation of the chiral condensate and the constituent quark mass in the Witten-Sakai-Sugimoto model are well known [19][20][21][22]50], and Fig. 3 underlines them. Nevertheless, one of our main points will be that the phase structure obtained in the present approach is very similar to the D3-D7 approach, where there is a straightforward concept of a chiral condensate and where the chiral condensate is zero in the chirally symmetric phase. This will be particularly obvious from the phase diagrams in the next subsection. In Fig. 4 we choose a fixed value of the mass parameter and plot the quark number density as a function of the chemical potential for various different temperatures in the left panel, and show the behavior of the (logarithm of the) pressure P = −Ω for one selected temperature in the right panel. The density curves can be understood as follows: at zero temperature (rightmost curve) the LTQ phase exists for arbitrarily large chemical potentials. As we increase the temperature, a point where the LTQ phase connects with the HTQ phase (red and blue curves in the figure), appears at largeμ and moves all the way down toμ = 0. In doing so, it bends and drags the LTQ curve backwards. Once the origin is reached, the LTQ curve detaches from the HTQ phase and connects at two points to the mesonic solution (n = 0), rather than interpolating between the mesonic and HTQ phases. For even larger temperatures (not shown here) this branch disappears and the LTQ solution ceases to exist. In the right panel we demonstrate that for certain temperatures there is a first-order LTQ-HTQ transition. For a better distinguishability of the curves, we have normalized all free energies by the free energy of the HTQ phase. We clearly see that there is a continuous transition from the mesonic phase to the LTQ phase, followed by a first-order transition to the HTQ phase.

Speed of sound
We have shown that there can be a continuous transition between the LTQ and HTQ phases. This continuity can be understood from the geometry in the bulk, but also manifests itself in the thermodynamic quantities considered so far. In particular, the quark number density is (within the very small numerical uncertainties) continuous across the transition, i.e., for all we know so far, the transition could be a smooth crossover. To analyze the smoothness of the transition, let us therefore compute second derivatives of the free energy. The speed of sound combines all possible second derivatives with respect to temperature and chemical potential and thus is a suitable quantity for our purpose. Comparable results of the speed of sound at nonzero temperature and chemical potential have been obtained in holography [31] as well as in NJL models [51,52]. The speed of sound c s is defined as where s is the entropy density and w = µn + st the enthalpy density, all partial derivatives with respect to t are taken at fixed µ and vice versa, and the derivative of the pressure P = −Ω with respect to the energy density = Ω + µn + st is taken at fixed entropy per particle s/n. In the definition (4.2) we have already replaced all thermodynamic quantities by their dimensionless counterparts. It is easy to check that this does not change the result for c s since all dimensionful constants and numerical prefactors cancel. We have computed the speed of sound using both expressions given in Eq. (4.2) and found that the numerically cleanest result is obtained by using the right-hand side with all second derivatives evaluated as far as possible in a semi-analytical way. These semi-analytical expressions (analytical up to a numerical integration) are derived in appendix A.
The numerical results are shown in Fig. 5. In contrast to the plots in the previous subsection, here we only show the results for the stable phases, i.e., metastable and unstable branches are omitted. Each curve thus contains contributions from two or three of the mesonic, LTQ, and HTQ phases, separated by phase transitions. We see that for large temperatures or chemical potentials the speed of sound approaches the massless result, as expected. This is a further check of our results since the speed of sound of the chirally symmetric phase in the massless limit can be computed from a completely different, much simpler expression [15] straight branes: c 2 s = 2 5 u T n 2 + u 5 T (n 2 + 5u 5 T ) + µn(n 2 + 6u 5 T ) where n and µ are related by In particular, all curves approach c 2 s = 2/5 for largeμ at fixedt and c 2 s = 1/6 for larget and fixedμ. Due to the inherent scale M KK and the lack of asymptotic freedom in this version of the model, these asymptotic values are different from the scale-invariant value c 2 s = 1/3 that is assumed asymptotically in QCD. Let us comment on the various phase transitions visible in Fig. 5. Thet = 0 curve in the left panel includes a mesonic-HTQ first-order transition which manifests itself in a discontinuity of the speed of sound. Thet = 0.15 curve shows a similar feature, but in this case it is a mesonic-LTQ transition with a smaller discontinuity. Since the numerics become challenging in the LTQ phase at very small n, it is difficult to compute the values of the speed of sound along the almost vertical segment of thet = 0.15 curve, i.e., the size of the jump is difficult to predict. More importantly, the same curve also shows a small, barely visible, discontinuity as we move from the LTQ to the HTQ geometry. In the vicinity of this transition, the numerics of the HTQ phase is completely under control. The reason is that, as mentioned above, the HTQ solution formally continues beyond the point x 4 (u T ) = 0 to unphysical, negative values of x 4 (u T ). Therefore, the point x 4 (u T ) = 0 is mathematically not very special. The numerics of the LTQ phase is somewhat more difficult as one approaches the point u c → u T , but the results for the speed of sound suggest no particular problems either. Therefore, we can be confident that this small discontinuity is not an artifact of the numerics. As a result, the LTQ-HTQ transition shown here is not a smooth crossover but a second-order phase transition, although with a very small discontinuity in the second derivatives of the thermodynamic potential. The same feature is visible in theμ = 0.8 curve in the right panel. Theμ = 0.3 curve in this panel shows the same mesonic-LTQ transition as in the left panel, but now followed by a first-order LTQ-HTQ transition, as discussed in the context of Fig. 4.

Phase diagrams
The main result of this paper is shown in Fig. 6. This figure shows four different phase diagrams, which summarize the various phase transitions discussed in the previous subsections. Before we proceed to the discussion of the results, let us briefly explain how we have calculated the phase transition lines. It would obviously be very tedious to determine the phase transitions point by point by calculating the free energies for all possible phases as functions ofμ ort. Moreover, it would then be easy to miss certain features of the phase structure, unless the free energies are computed on a very fine grid inμ andt, which would be very laborious. It is thus much more advantageous to compute the phase transition lines in a more direct way. For instance, for a first-order phase transition curve we have solved a system of equations simultaneously that combines the equations of the two respective phases that coexist at this phase transition line plus the condition that they have the same free energy. All curves in Fig. 6 have been obtained in this direct way, i.e., by suitably combining the stationarity equations and various conditions for the free energies.
The figure contains four phase diagrams for four different values of the (rescaled) mass parameterα, with the quark mass increasing from top to bottom rows. As in the previous section,λ = 15. The four phase diagrams in the left panels are duplicated in the right panels, but with additional auxiliary curves that are helpful for the understanding (and for reproducing) the phase structure since they indicate the regions of unstable and metastable solutions. In particular, these curves define the regions where each of the three phases exist and keep track of the rich behavior of the LTQ solution: • The HTQ phase exists everywhere above the (red) curve labeled by r. On this curve, x 4 (u T ) = 0 and it connects to the LTQ phase. This solution is always single-valued.
The four values of the mass parameter are chosen to represent prototypical examples of the four different topologies of the phase diagram we have found. For small current quark masses the phase diagram is qualitatively the same as in the chiral limit: there is a discontinuous transition between the mesonic phase, which has zero quark number, n = 0, to the phases with nonzero quark number n > 0. In contrast to the chiral limit, the n > 0 region now contains the LTQ phase, whose area grows as we increase the current quark mass. We have depicted the transition between LTQ and HTQ phases by a dotted line. This indicates that the transition is continuous, i.e., density and entropy do not jump when we cross from one phase to the other. However, we have seen in the calculation of the speed of sound that second derivatives have a (very small) discontinuity. This discontinuity might be a consequence of constructing the LTQ phase with the help of string sources, which possibly is an approximation to a phase with a more complicated flavor brane embedding, as already discussed at the end of Sec. 3.2. We cannot exclude that a more complete treatment removes this second-order discontinuity and results in a smooth crossover.
As we move to larger values of the current quark mass, a first-order phase transition line develops at intermediate temperatures. The endpoint of this line is somewhat difficult to calculate because the free energies of the phases are extremely close to each other in the vicinity of this point. Therefore, it is also difficult to calculate the critical value of the mass parameter where this first-order line starts to appear. Our numerics suggestα 0.007. Increasing the mass parameter further, the original first-order chiral transition line breaks up, and a continuous transition between the mesonic and the LTQ phases appears. The exact location of the resulting endpoint of the first-order mesonic-LTQ transition is again   Figure 7. Different topologies of the phase diagrams in the plane of rescaled 't Hooft coupling λ and quark mass parameter, classified according to the mesonic-LTQ transition, which can be discontinuous everywhere ("1st"), discontinuous at small temperatures and continuous at large temperatures ("1st & 2nd") or continuous everywhere ("2nd"). The shaded area is excluded since our approximation breaks down. The 4 (red) dots mark the parameter pairs chosen for the phase diagrams in Fig. 6. The translation from the mass parameterα to the more physical quantity f 2 π m 2 π is obtained with the help of Table 1.
prone to some numerical uncertainties. Eventually, at even larger mass, the mesonic-LTQ transition becomes continuous throughout.
The "evolution" of the phase diagram discussed so far is obtained at fixed (rescaled) 't Hooft couplingλ, and one may ask how the phase structure changes if alsoλ is varied. (Keeping in mind that our classical gravity approximation is only valid at strong 't Hooft coupling λ and results for small values of the coupling are an extrapolation.) The result is shown in Fig. 7, which indicates the different topologies of the phase diagram in thẽ λ-α-plane. Recall that since we were able to absorb the other model parameters and M KK in our rescaled variables,λ andα are the only independent parameters, and thus Fig. 7 presents a complete systematic study of the phase structure. This figure contains an unphysical region for large values of the mass parameter, shown as a shaded area. As we enter this area, a region in thet-μ plane appears where there is no stable phase in our approximation. The reason is that the ∪-shaped solution ceases to exist for small temperatures (for all temperatures as we move deeper into the shaded area), as already pointed out in Sec. 4.1. This shortcoming of our approach is no surprise since we have only added a term linear in the current quark mass to the action.
For the horizontal axis in this figure we have translated the mass parameterα into the more physical quantity f 2 π m 2 π /L −4 via Table 1 with N c = 3. For a rough estimate, let us use L 10 −3 MeV −1 , which was obtained from a fit to nuclear matter properties at saturation density [15]. Then, with the actual QCD vacuum values f π 93 MeV and m π 140 MeV we have f 2 π m 2 π /L −4 1.7 × 10 −4 . Atλ = 15 this corresponds toα 8 × 10 −4 , slightly smaller than even the smallest mass parameter used in Fig. 6. Thus all interesting features of our phase diagram occur in the heavy QCD regime. We may also use this fit to get a rough idea of the scales in the phase diagrams. With the help of Table 1 we find that for the chiral limit the µ = 0 phase transition is at about T c 150 MeV, while the t = 0 (discontinuous) phase transition occurs at a quark chemical potential of about µ c 520 MeV. The heaviest quark mass parameter considered in Fig. 6 corresponds to f 2 π m 2 π being about 250 times its QCD value. In this case the corresponding critical values are T c 240 MeV and (now for a continuous transition) µ c 1200 MeV.
It is tempting to compare our results to lattice QCD, where at large quark masses and large coupling the phase structure can be computed even for nonzero chemical potentials [10,32]. If we take our heaviest case and ignore the fact that our continuous phase transitions have higher-order discontinuities, the phase diagram has a intriguing similarity to the heavy QCD diagram, see for instance Fig. 1 of Ref. [10]. The similarity becomes even more striking when we go to our intermediate regime where there is a low-temperature discontinuous transition with a critical endpoint. However, this comparison has to be taken with some care for at least two reasons. Firstly, our calculation does not include baryons. Therefore, our low-temperature mesonic-LTQ transition that ends in a critical point is different from the liquid-gas transition in Fig. 1 of Ref. [10]. Secondly, heavy QCD becomes pure Yang-Mills in the limit of infinitely heavy quarks. In this limit, the first-order transition at large temperature is a strict deconfinement transition. In our approach, however, we have decoupled the gluon dynamics such that our high-temperature transition does not seem to know anything about confinement. Nevertheless, the resemblance of the phase diagrams with heavy QCD is striking and suggests further studies for a better understanding of the parallels and differences to QCD, most notably by including baryons into our holographic calculation.
A more concrete comparison can be made to the phase diagrams obtained in holographic studies with a D3-D7 setup. We have already pointed out the very close analogy of our three flavor brane configurations to the corresponding D3-D7 geometries. This analogy is also borne out in the phase diagrams, as one can see for instance by comparing our results to Fig. 5 in Ref. [26] (where a nonzero magnetic field was introduced in order to break conformal symmetry). This particular figure is obtained with zero current quark mass, but in that model the quark mass does not seem to change the topology of the phase diagram except for turning the chiral phase transition into a smooth crossover. The first-order transition, including its critical endpoint, as well as the overall structure is qualitatively the same as in our heaviest case, see lowest panels in Fig. 6, if the meson-melted D3-D7 phase is identified with our LTQ phase. This very close resemblance is somewhat surprising since, as already discussed, in the present model the identification of an order parameter for chiral symmetry breaking is much less straightforward than in the D3-D7 approach. Turning this argument around, we can take the result as a reassurance for our approach, despite the open conceptual questions regarding the implementation of the chiral condensate in the Witten-Sakai-Sugimoto model.

Summary and outlook
We have included a nonzero current quark mass in the study of the phase diagram in holographic QCD. Working in the deconfined geometry and the decompactified limit of the Witten-Sakai-Sugimoto model, we have added a mass correction to the action analogous to chiral perturbation theory. To this end, we have made use of earlier works that identified the chiral condensate with the expectation value of a non-local, gauge-invariant open Wilson line operator. The holographic dual of this expectation value corresponds to a worldsheet instanton, and we have included this contribution for the first time in a fully consistent way, calculating its effect on the embedding of the flavor branes in the background geometry. We have identified three different embeddings. Most notably, we have pointed out that the current quark mass stabilizes a brane configuration with string sources, which previously had only been discussed in the chiral limit, where it is never stable. Evaluating the classical equations of motion and the stationarity conditions of the on-shell action, we have systematically studied the equilibrium phase structure in the plane of temperature and baryon chemical potential. The effect of heavy quarks is to break up the chiral phase transition line and introduce a critical endpoint at high temperatures. At intermediate values of the quark mass there is an additional endpoint at lower temperatures. This phase structure is qualitatively similar to heavy QCD on the lattice. We have also pointed out that the phase structure as well as the geometric structure of the different embeddings is strikingly similar to a different holographic approach based on the D3-D7 setup. This similarity is much less obvious in the chiral limit.
The main significance of our work is to make the study of equilibrium phases in the Witten-Sakai-Sugimoto model more realistic. Therefore, our results can be used as the basis for further studies of strongly-coupled hot and dense matter. The most obvious extension is to include baryonic matter, possibly as a first step from simple approximations based on pointlike baryons [8], or from more sophisticated approximations using nonzerowidth, interacting instantons [15]. This extension would be very interesting for the phase structure itself but also for more specific questions such as the quark-hadron continuity at low temperatures [15] and for a realistic equation of state for dense matter, to be applied to the physics of neutron stars.
Our work also suggests several other directions for future research. For instance, we have only considered the case of a degenerate quark mass for all flavors, such that the number of flavors played no important role. The generalization to different quark masses is straightforward and could have non-trivial effects on the phase structure of the system. We have also pointed out and reinforced several conceptual problems related to the introduction of a current quark mass in the Witten-Sakai-Sugimoto model. This concerns for instance the identification of chiral condensate, constituent quark mass, and order parameter for chiral symmetry breaking. Our results for the phase diagram, which seem very sensible in comparison to other model approaches and even compared to actual QCD, suggest that one might simply go ahead with the present approach. Nevertheless, some conceptual progress would be desirable as well. It might be useful to repeat our study employing another approach that is suggested in the literature for the introduction of a mass term, based on a tachyonic scalar field. An improvement to our study on the conceptual level would also be to find a smooth flavor brane configuration that replaces our construction with string sources and a cusp in the embedding of the branes. derivative of Ω with respect to u c vanishes at the stationary point. Computing derivatives in this way, i.e., via L 0 , which is a functional of x 4 (u) andâ 0 (u), is equivalent to using the explicit form of the free energies in Sec. 3. Since x 4 (u c ), x 4 (∞) andâ 0 (∞) are held fixed and thus their derivatives with respect to u T vanish, and using For small temperatures, since u T ∝ t 2 , the entropy grows linearly with t. One easily finds that for the mesonic phase the entropy density s ∪ has the same form with n = 0. In this phase the speed of sound is simply given by For the HTQ phase the analogous calculation gives Next we compute the number susceptibility. For the LTQ phase let us define We compute ∂µ ∂n at fixed t along the surface given by the constraints f A = f k = 0. These constraints are the first two relations in Eq. (3.22). To this end, we first write on account of the implicit function theorem (A.9) The chemical potential µ is a function of n, u c , and A via the third relation in Eq. (3.22 (Temporarily, we have employed the notation of a total derivative on the left-hand side, to avoid confusion with the explicit derivative on the right-hand side.) All partial derivatives can be taken semi-analytically by taking the derivatives of the integrands in the various integrals. The result is a complicated combination of numerical integrals which however is easily handled numerically. The number susceptibility ∂n ∂µ is then simply obtained by the inverse of Eq. (A.10).
It remains to compute the number susceptibility in the HTQ phase. This is somewhat simpler because here u c = u T and thus there is one fewer dynamical variable. We define where the explicit derivative of µ with respect to n is obtained from the second relation in Eq. (3.28). This concludes the collection of all derivatives for all three phases needed for the calculation of the speed of sound according to the strategy laid out at the beginning of this appendix.