Holographic QCD with Dynamical Flavors

Gravity solutions describing the Witten-Sakai-Sugimoto model of holographic QCD with dynamical flavors are presented. The field theory is studied in the Veneziano limit, at first order in the ratio of the number of flavors and colors. The gravity solutions are analytic and dual to the field theory either in the confined, low temperature phase or in the deconfined, high temperature phase with small baryonic charge density. The phase diagram and the flavor contributions to vacuum (e.g. string tension and hadron masses) and thermodynamical properties of the dual field theory are then deduced. The phase diagram of the model at finite temperature and imaginary chemical potential, as well as that of the unflavored theory at finite theta angle are also discussed in turn, showing qualitative similarities with recent lattice studies. Interesting degrees of freedom in each phase are discussed. Covariant counterterms for the Witten-Sakai-Sugimoto model are provided both in the probe approximation and in the backreacted case, allowing for a standard holographic renormalization of the theory.


Introduction
The Witten-Sakai-Sugimoto (WSS) model [1,2] is the most studied top-down example of "holographic QCD". Its realization in string theory describes, in the low energy limit, a four dimensional SU (N c ) non-supersymmetric Yang-Mills theory coupled to adjoint massive matter and N f chiral massless fermions -quarks in the fundamental representation. Its popularity is due to its ability in describing in a very simple and calculable way many physical properties interesting for QCD, such as confinement, chiral symmetry breaking, confinementdeconfinement transition, etc. 1 In fact, the WSS model is the holographic theory closest to (planar) QCD, 2 although it is not its precise dual.
In the original form of the Witten-Sakai-Sugimoto model the quarks in the fundamental representation are quenched, non-dynamical degrees of freedom. In this paper we construct holographic solutions corresponding to dynamical fundamental matter, at first order in the Veneziano limit. That is, the gravitational solutions we present provide an accurate description of the field theory in the limit of large number of colors and flavors, N c 1, N f 1, with fixed ratio f ∼ N f /N c , at first order in the small f -expansion. We study the model both in the confined and deconfined phase (below and above a critical temperature) and we also include a finite, small baryonic charge density. Our solutions describe a configuration where the flavor symmetry is completely Abelian -U (1) The gravity solutions allow for the study of effects of dynamical flavors on the observables of the theory. As illustrative examples, we first briefly consider the string tension and the mass spectrum of some hadrons in the confined case. Then, we re-derive the thermodynamics of the solutions (previously studied in the quenched limit in [3,4,5]) and discuss the phase diagrams (in a section hopefully accessible to the non expert on holography). Our aim is to compare the holographic results to those obtained in lattice QCD. We first consider the flavor dependence of the deconfinement temperature 3 and its behavior as a function of the baryon chemical potential, obtaining qualitative matching with lattice QCD results. Then we focus on the phase diagram of the model at finite temperature and imaginary chemical potential, as well as on that of the unflavored theory at finite θ angle. The holographic model is observed to provide a natural realization of a qualitative "duality" between these two phase diagrams, recently observed in lattice (quenched) QCD [6]. In the latter paper, in fact, the phase diagram at finite temperature and θ angle is observed to have the same qualitative features as an "inverted" (Roberge-Weiss) phase diagram at finite temperature and imaginary chemical potential. In the holographic model this is mainly due to a symmetry of the background under the exchange of two circles in the geometry. The holographic model also allows for an immediate proposal for the relevant degrees of freedom in the two phases: instantons in the deconfined case and Euclidean baryon vertices in the confined case.
Coming to some technical details, we construct the solutions by making use of the smearing technique [7,8] (see also [9] and [10,11] for reviews). The WSS model consists in adding "flavor" D8-branes to the background generated by D4-branes wrapped on a circle [1]. This background is dual, in the low energy limit, to a Yang-Mills theory with adjoint massive matter. The D8-branes, which describe the flavor sector, are placed at a fixed value of the wrapped circle mentioned above [2]. Studying their backreaction on the background in the localized case is a daunting task and it has been performed only in a particular limit in [12].
We are able to take into account the backreaction of the D8-branes, and so to consider the dynamics of fundamental matter, by homogeneously smearing a large number of D8-branes on the wrapped circle, 4 such that the isometries of the original background are preserved. This configuration greatly simplifies the setting, reducing the equations to be solved to ODEs instead of PDEs. The system is still extremely complicated: it is non-supersymmetric -giving a set of coupled second order equations -and the finite charge density implies the presence of highly non-linear terms. Nevertheless, surprisingly enough, the system admits analytic solutions in the limit of small flavor backreaction and small charge density. The parameters of the solutions can be chosen in such a way that the backgrounds are completely regular in the dual IR regime. In the dual UV regime, on the contrary, they present a non-removable divergence, due to the presence of a Landau pole in field theory. 5 The latter is reflected in the holographic running coupling too.
Last but not least, in order to compute the free energy of the system, we holographically renormalize the theory. As far as we know, despite a huge number of studies of the WSS model, its holographic renormalization is not known in a covariant way so far. Thus, we provide covariant counterterms for the WSS model, filling the gap in the literature. The form of the counterterms is the same both in the probe approximation and in the backreacted case (at leading order in f ).
The paper is organized as follows. In the next section we provide a very short review of the WSS model in the quenched case. Then, in section 3 we construct the gravitational solution with backreacting D8-branes in the dual confined phase. Some physical quantities, such as the running coupling, string tension and the baryon vertex and vector meson mass spectra, are briefly discussed as well. The holographic renormalization of the WSS model is described in section 4. In section 5 we present the gravitational solution dual to the deconfined phase with finite charge density, and we derive its thermodynamics in section 6. Section 7, which is written also for non experts in holography, presents the discussion of the critical deconfinement temperature , the phase diagrams and some interesting degrees of freedom of the theory. We conclude in section 8 with some comments on the solution, its limitations and potential applications. In the appendix we write the solution in the particular limit considered in [13].

A review of the Witten-Sakai-Sugimoto model
In [1], a large N c non-supersymmetric Yang-Mills theory in 3 + 1 flat space-time dimensions was obtained by Witten as the low energy limit of a Kaluza-Klein (KK) reduction of a particular 4 + 1 dimensional SU (N c ) Yang-Mills theory coupled to massless adjoint scalar and fermionic matter fields. This theory describes the low energy dynamics of open strings whose end-points are attached (with Dirichlet boundary conditions) to a stack of N c parallel D4-branes (4 + 1-dimensional hyperplanes) placed in an ambient ten dimensional Minkowski space-time. The reduction is realized by compactifying the theory on a circle S x 4 of length β 4 along the fourth space direction x 4 and choosing periodic (resp. anti-periodic) boundary conditions for bosons (resp. fermions). In this way the massless modes at energy E 1/β 4 are those of a 3+1 dimensional SU (N c ) Yang-Mills theory. The other modes get masses M KK of the order of 1/β 4 . If we denote by T s the string tension of the confining 3 + 1 dimensional Yang-Mills theory, the low energy theory can be decoupled from the Kaluza-Klein modes provided T s /M 2 KK ≡ 2λ 4 /27π 1. Here, λ 4 is the 4d UV 't Hooft coupling. Unfortunately, Witten's model in the most interesting λ 4 ∼ 1 regime has no known simple description. However, we can obtain many detailed informations in the λ 4 1 regime, where the theory (in the large N c limit) is conjectured to have a holographic dual description in terms of a classical theory of gravity on a background which arises as the near-horizon limit of the one sourced by the N c D4-branes. This background has the topology of a product Here R 1,3 is the 3 + 1 dimensional flat Minkowski space-time, while R u denotes the real semi-axis spanned by a radial coordinate u which is roughly the geometric counterpart of the Renormalization Group energy scale in the dual field theory. The (u, x 4 ) subspace is a cigar, with the S x 4 circle smoothly shrinking to zero size at a finite value u 0 of the radial coordinate. 6 Finally, S 4 is a compact four-dimensional sphere, whose isometry group SO(5) is holographically mapped into a global symmetry group under which the massive Kaluza-Klein fields (which are not decoupled in the large λ 4 regime) rotate. Despite being in a regime where the interesting Yang-Mills part is intrinsically coupled with spurious Kaluza-Klein modes, the λ 4 1 theory displays confinement and (once coupled to chiral massless quarks) chiral symmetry breaking at T = 0.
The introduction of N f chiral fundamental massless quarks in Witten's model can be achieved by adding suitably embedded N f D8 −D8-branes (which are localized in the x 4 Figure 1: The cylinder and cigar sections of Witten-Sakai-Sugimoto's gravity solution at T > T c (deconfined phase) and T < T c (confined phase). Here u is the radial coordinate (the geometric counterpart of the RG energy scale in the dual field theory) and x 0 and x 4 are Euclidean time and space directions compactified on circles of length 1/T and 2π/M KK respectively. The N f D8-branes (resp.D8-branes) account for the presence of N f right-handed (resp. left-handed) massless quarks in the dual quantum field theory. These branes connect in a "U"-shaped configuration at T < T c : this provides a simple geometrical realization of chiral symmetry breaking. direction) in the background, as described by Sakai and Sugimoto in [2]. In fact, quark fields arise in the massless spectrum of the open strings stretching between the "color" D4-branes and the "flavor" D8-branes. Neglecting the flavor brane backreaction on the background amounts on treating the quarks in the quenched approximation. The flavor branes support a U (N f ) × U (N f ) gauge group which is holographically mapped into the global flavor symmetry of the dual field theory. Chiral symmetry breaking is simply realized in the model: at T = 0 the flavor branes will find energetically convenient to connect into a "U"-shaped configuration, see figure 1. This will automatically break the original flavor symmetry group into a diagonal SU (N f ) × U (1) B subgroup, where U (1) B is the baryon symmetry. The axial symmetry U (1) A , instead, is broken by the anomaly.
Going to finite temperature amounts on compactifying the (Euclideanized) time direction x 0 on a circle S x 0 of length β = 1/T . There are two allowed gravity backgrounds where this is realized. A first one has the same structure as that at T = 0. The (x 0 , u) subspace is a cylinder while the (x 4 , u) one is still a cigar. The full background is given by Here µ = 0, 1, 2, 3. R measures the curvature radius of the background in string units, R = (πg s N c ) 1/3 l s , where g s is the string coupling and l s is the string length. Moreover dΩ 2 4 is the metric of the transverse S 4 , φ is the dilaton, F 4 is a Ramond-Ramond four-form and ω 4 is the volume form of S 4 . The 't Hooft coupling of the 4d theory in the UV is given by The shrinking of the S x 4 circle at u = u 0 is evident from the fact that f (u 0 ) = 0. The (u, x 4 ) subspace is a cigar, i.e. there are no conical singularities at u = u 0 provided the relation 9β 2 holds. The x 4 periodicity β 4 is related to the mass scale M KK of the Kaluza-Klein modes by . A second allowed solution (which thus has the same u → ∞ asymptotic as the one above) is instead a black-hole background 7 with Euclidean metric given by Here a = 1, 2, 3. The (x 0 , u) subspace is a cigar -with the S x 0 circle smoothly shrinking to zero size at u = u T (the position of the horizon) -provided we identify The (x 4 , u) space is instead a cylinder. Notice that the two allowed background metrics are related by an exchange of S x 0 with S x 4 and of u T ∼ T 2 with u 0 ∼ M 2 KK (notice in turn that these relations suggest a radius/energy relation of the form u ∼ E 2 ).
The background which minimizes the Euclidean on-shell gravity action S will provide the dominant contribution to the gravity partition function in the classical limit Z ∼ e −S . It turns out that the black hole solution is the energetically preferred one when T > T c , where the critical temperature T c,0 is given by (2.5) 7 It has been claimed in [14] that the black hole solution is not smoothly connected to the deconfined phase of pure YM in 4d for λ 4 1. The "correct" background in this sense is not explicitly known. We thus study the original construction as a model for the confinement/deconfinement transition.
At lower temperatures the other solution turns out to be energetically preferred.
Each background solution is dual to a different phase of the corresponding field theory. At T > T c the theory turns out to be in a deconfined regime, while at T < T c the theory is confining. This can be explicitly checked by e.g. a holographic computation of the Polyakov loop.
The holographic correspondence maps the difference between the on-shell Euclidean gravity actions on the two backgrounds with the free energy difference, among the two phases, in the dual large N c , large λ 4 quantum field theory: ∆S = β∆F . One can thus compute ∆F and study the phase diagram in a relatively simple way. In particular one can immediately verify that there is a first order phase transition at T = T c .
When the D8 −D8-branes are placed at antipodal points on S x 4 (the configuration we consider in this work), they fall into the horizon of the black hole in the T > T c phase, see figure 1. Thus, they do not reconnect anymore and the flavor symmetry group remains U (N f ) × U (N f ): chiral symmetry is restored in the deconfined case and the transition temperature T c coincides with the confinement/deconfinement one [3]. When the D8−D8-branes are not antipodal, deconfinement and chiral symmetry restoration generically happen at different critical temperatures [3]. Though this is an interesting possibility to consider, we leave the analysis of the corresponding setup for the future.

The backreacted solution in the confined phase
In [2] the flavor D8-branes were treated in the probe approximation. Here we want to take into account their backreaction on the background to first order in N f /N c . The trick we use in order to avoid having to deal with partial differential equations coupled to delta-function sources, it to consider a setup where a large number N f of D8-branes are homogeneously smeared along the transverse x 4 circle [7,8,10,11]. In this section we consider the whole model at T < T c and at finite quark chemical potential for the flavor fields.
In our analysis (both in the confined and in the deconfined case) we will keep the length of the x 4 and x 0 circles fixed. The 't Hooft coupling λ 4 , "geometrically" defined from the reduction of the 5d UV completion of the model on x 4 , will thus be held fixed too.

Action and ansatz
Our string frame metric ansatz for the T = 0 case is the same one used in [3] 8 ds 2 = e 2λ (−dt 2 + dx a dx a ) + e 2λ dx 2 4 + l 2 s e −2ϕ dρ 2 + l 2 s e 2ν dΩ 2 4 , where λ,λ, ϕ, ν are functions of the radial coordinate ρ and the x 4 -coordinate is compactified on a circle of length When T = 0 what we need to do is just to compactify the Euclidean time on a circle of length β = 1/T . Finally, the function ϕ is related to the dilaton φ through the defining equation The relevant action from which the equations of motion follow is 9 The last part of this action arises as the DBI contribution of N f U-shaped D8-branes which are homogeneously smeared on the transverse x 4 circle (whose metric component is g 44 and whose length is 2π/M KK ) and reach the tip of the cigar. 10 At leading order in N f /N c further contributions from the bulk gravity fields will not enter the analysis and will hence be neglected (they are zero on-shell on the ansatz we use). 11 The smeared DBI action is put on-shell w.r.t. the embedding coordinate x 4 = x 4 (ρ), which, in the confining case and for the antipodal configuration we are considering, satisfies the equation of motionẋ 4 = 0. In order to account for the presence of two branches at two antipodal points on the x 4 circle, we count the integration over the radial coordinate two times. Moreover, we put on-shell also the U (1) gauge field on the branes, which is the holographic dual of the U (1) B current. A non trivial value of the electric component of the U (1) field strength F is related to a finite baryon density configuration. In the present case F ≡ dA = 0: this in fact the relevant solution in the confining phase [4] at small baryon chemical potential µ and if no explicit sources are introduced. 12 This corresponds to field theory configurations, which are encountered also in the QCD phase diagram, with finite quark chemical potential and zero baryon density.
In the probe approximation, with localized flavor branes, there are other possible relevant configurations of the gauge field. For example, in the deconfined case the preferred configuration above a critical value of the charge density is spatially modulated [18]. In the confined phase one can consider non-Abelian configurations and add to the setup D4-branes wrapped on S 4 which are instantons on the D8 worldvolume and act as sources (baryon vertices) of 9 Notice that in the following, since the g s factor is already contained in the on-shell value of e φ as in the original Witten-Sakai-Sugimoto model, the D8-brane tension T 8 will be defined as T 8 = (2π) −8 α −9/2 . 10 Notice that we smear the U-shaped configuration along of the whole x 4 circle. 11 See e.g. [15,16,17] for the similar situation in the D3 − D7 system. 12 The solution corresponds to a constant temporal component of the brane gauge field, A t = µ, which decouples from the remaining equations of motion. Notice that a constant value for the temporal component of the gauge field is allowed since, in the confining phase, the temporal circle does not shrink. the U (1) B field (see for example [5]). We will not consider these involved configurations in our context. Implementation of the ansatz above (plus the usual ansatz in (2.1) for F 4 ) gives the following 1d action which has to be supported by the zero-energy constraint as the constant arising from the quantized F 4 flux through S 4 and as the one related to the number N f of flavor branes. Finally the overall volume factor reads where V 3 is the infinite 3d-space volume and V S 4 = 8π 2 /3.

Equations of motion
The equations of motion following from the previous action, which we re-arrange so that the dilaton φ appears instead of ϕ, arë In what follows it will be useful to define a new radial coordinate where u 0 is the minimal value of the original Witten-Sakai-Sugimoto radial variable u, i.e. the position of the tip of the (x 4 , u) cigar in the unflavored case.
We are going to look for a perturbative solution of the above equations at first order in the effective parameter which weighs the flavor contribution to the action. In particular we will expand all the functions in the form is our expansion parameter and Hooft coupling at the UV scale set by M KK . The zero-th order unflavored solutions in these coordinates can be read from [3] and are given by In the unflavored case, the variable r is related to the standard Witten-Sakai-Sugimoto radial variable by The equations of motion for the first order terms read (derivatives are with respect to r now) From these we see thatλ and A 1,2 , B 1,2 are integration constants.

Solution and asymptotics
Plugging the expressions (3.18, 3.19) in the remaining equations and integrating, we find the solution given in terms of (generalized) hypergeometric functions. All in all we have eight integration constants A 1,2 , B 1,2 , C 1,2 , M 1,2 . Some of them are fixed by physical requirements. The zero energy constraint (3.6) is satisfied, to first order in f , provided the condition holds.
Other constraints arise by requiring regularity at the tip of the (x 4 , r) cigar, which corresponds to the r → ∞ limit. In this case the corresponding (IR) behavior of the various functions above can be easily extracted working with the variable x = e −3r/2 , which goes to zero in the limit. As a result we find the following IR asymptotics Regularity at the tip of the cigar can be achieved imposing the following conditions which can be further restricted by requiring in order to drop all the IR logarithmically divergent factors from the metric components and the dilaton; this is not strictly necessary: the logarithmically divergent factor inλ 1 is anyway subleading, in the small f limit, w.r.t. the analogous term inλ 0 . It is notable that if condition (3.24) is satisfied, the constraint (3.22) following from the zero energy condition is automatically satisfied too. The UV (r → 0, i.e. x → 1) behavior of the above functions is given by It is interesting to notice that the subleading UV divergences, that is the terms diverging as (1 − x) −1/6 , do not depend on the integration constants. We interpret these as the dual of the "universal" terms related to the addition of the flavors to the original unflavored action. As in the D3 − D7 case [9], the flavor source term induces a Landau Pole in the UV. This is the reason why these subleading divergences cannot be killed by a choice of integration constants.
The combinations of integration constants appearing in the UV asymptotics, can instead be interpreted as corresponding to sources or VEVs of other gauge invariant operators. The sources can be turned off by a suitable choice of the integration constants. We defer the field/operator analysis for the future, so we do not know yet what are the combinations of λ 1 ,λ 1 , φ 1 , ν 1 corresponding to gauge invariant operators. Nevertheless, we are interested in switching off at least the most divergent terms in (3.26), hence we impose the prudent condition A supplementary option is to make sure that all the subleading terms are turned off, which requires the condition These conditions, together with (3.29), imply that In the following we will consider the family of solutions for which the IR regularity conditions (3.24), and the conditions (3.29) killing the leading UV divergences hold. On this family we thus have We will take A 1 , A 2 , C 2 generic, having in mind their special values (3.25), (3.31) along the analysis.

Some physical properties
In this section we give a flavor of the physical effects which can be studied with the backreacted solution. Let us begin by analyzing how the flavor contributions modify the relation between the radial parameter u 0 and the mass M KK . This relation comes from the requirement that the (x 4 , r) cigar closes smoothly at the tip (which is at r → ∞) and reads For the special choices (3.25), (3.31),

The running coupling and the UV Landau Pole
The running gauge coupling can be obtained by examining the action of a probe D4-brane wrapped on the x 4 circle (see e.g. [19]). It reads (using the (3.34) hence, differently from the unflavored case (where it goes to the constant "geometrical" value g 2 Y M ), the running coupling tends to diverge in the UV, signaling the presence of a Landau pole. The gauge theory analysis can thus be safely performed only in the IR, that is at scales Note that in the perturbative analysis we are considering (where f 1) the UV Landau Pole essentially coincides with the asymptotic value x → 1.

The string tension
Using (3.33), the string tension turns out to be given by Notice that this does not depend on C 2 . For the special choice of integration constants (3.31), 3A1 − A2 − 28k ≈ 1.13. In this case one could conclude, naively, that the effect of the dynamical flavors is that of increasing the string tension. However, as we will also discuss in the following, care as to be taken when comparing theories with different number of flavors.

Baryon mass
In the model at hand, a baryon vertex is identified [20] with a Euclidean D4-brane wrapped on S 4 and localized at the radial position corresponding to the deep IR of the dual field theory. Using the backreacted metric found above we can easily study how the dynamical flavors affect the mass of the baryon. The wrapped D4-brane action reads is the D4-brane tension and x = e −3r/2 → 0 is the IR value of the radial variable. Using our solution we thus get that the baryon mass is given by For the special choices (3.31),

Vector meson mass spectrum
Following the same reasonings as in Section 3.3 of [2] (to which we refer the interested reader for details), the spectrum of the massive vector mesons on the backreacted background can be found by considering the fluctuation of a gauge field on a probe D8-brane. The massive vector mesons have masses given by where γ n are given by the equation 42) and the functions ψ n satisfy the normalization condition Normalizability forces us to choose an asymptotic condition such that the modes ψ n vanish in the UV. As we have already noticed, the flavored theory we are considering has a Landau Pole in the UV and hence the same condition has to be imposed at a UV cutoff. In the deep IR (i.e. at x = 0 or r → ∞) we impose the regularity conditions ψ = 0 or ψ = 0 which distinguish among even and odd functions. Correspondingly we have vector and axial-vector mesons. Modes with n odd (resp. even) correspond to the former (resp. the latter). As in [2] we have that the lightest mode γ 1 corresponds to a vector meson (the "ρ") with C = P = −1, the second one γ 2 corresponds to an axial-vector meson with C = P = +1 (the "a 1 (1260)") and so on. We solve equation (3.41) with the standard shooting technique. We do not perform a scan of the results as the parameters of the solution are varied, but consider one relevant case. Namely, the integration constants are fixed such that: i) all the IR logarithmical divergences are canceled (i.e. we enforce the special condition C 2 = 0); ii) the free energy coincides with the one obtained in the probe calculation (as expected, see the next section, formula (4.77)); iii) the source for the operator dual to the field λ 1 −λ 1 is turned off. The latter combination enters naturally the reduction of the ten dimensional action to five dimensions (where the unflavored metric is asymptotically AdS in the dual frame [21]), giving a well defined scalar in the reduced theory, so it is very likely to have a dual operator without further mixing with other modes. Its source depends only on the constant A 1 . Conditions ii) and iii) turn out to coincide with (3.31), giving The results for the vector meson spectra are as follows. The ratio of the masses of the first two ρ mesons is reduced by the flavor contribution, 13 giving a result closer to the phenomenological one than the unflavored theory. Namely, while experimentally m 2 ρ(1450) /m 2 ρ ∼ 3.5 and in the Witten-Sakai-Sugimoto model this ratio is around 4.3, in the flavored case with f = 0.02 (taken as a representative value) it is 3.7. Also the ratio of the masses of the first axial vector with the ρ meson is reduced by the flavor contribution, giving a result which is now more distant from the phenomenological one than the unflavored theory. That is, from experiments m 2 a 1 (1260) /m 2 ρ ∼ 2.51, in the Witten-Sakai-Sugimoto model this value is 2.39 and in the flavored case with f = 0.02 it is 2.37. As in the unflavored theory, the masses of higher mesons are quite larger than the experimental ones.
In the unflavored case, it is possible to look for approximate solutions of the equations (3.41) also using the WKB method as reported e.g. in [23]. The result should hold in principle only for higher modes, but it is "massaged" in order to account for the small n behavior too. It reads to be compared with the numerical results 0.67 , 1.6 , 2.9 , 4.5 obtained in [2]. In the flavored case the WKB analysis is harder. One can imagine that (especially for low lying modes) the meson mass ratios do not depend on the f parameter (as the f dependence should be encoded in the overall scale). If this is the case one could guess that where ξ is the flavored version of ξ 0 where I is a divergent integral. This indicates that a careful treatment of the solution is needed by introducing a cutoff and expanding the solution below this cutoff. We hope to come back to this issue in the future.

Holographic renormalization of the Witten-Sakai-Sugimoto model
In this section we present the holographic renormalization of the Witten-Sakai-Sugimoto model, both in the probe approximation and in the backreacted case. As far as we know, the covariant counterterms needed to renormalize the action of this very well known model are not present in the literature. In the paper [24], where the WSS model in an external magnetic field is considered, counterterms for the probe D8-branes are provided, but they are not written in a covariant way. A contingent reason why we need to holographically renormalize the theory is that we want to compute the free energy of the model. This is related to the renormalized on-shell Euclidean gravity action through the holographic formula F = T S ren E,on−shell . Here where S E is the Euclidean version of the bulk+smeared-D8-brane action Notice that we work in string frame. The standard Gibbons-Hawking term reads 14 where h is the determinant of the boundary metric (in our coordinates it is the slice of the 10d metric at fixed r = with → 0 being the UV boundary), K is the trace of the extrinsic curvature of the boundary, which in our case is explicitly given by and is an outward pointing unit normal vector to the boundary at r = (the minus sign originates from the fact that our coordinate r decreases towards zero when approaching the boundary).
The "bulk" counterterm is a volume boundary term, already introduced in [25] S bulk c.t. = 1 Evaluating the above terms on our solution we get In the unflavored case, the bulk counterterm precisely cancels the 1/ divergence from S E + S GH so that consistently with equation (E.13) in [25]. Correspondingly the free energy density reads (using (2.2)) 15 The novelty in our analysis at this point is that in the flavored case there are (apart from novel constant terms) two classes of novel divergences (in −1/6 and −7/6 ) and the counterterm we have introduced above is not enough to cancel them We thus need to introduce new counterterms S D8 c.t. for the D8-brane contribution.

The probe approximation
In the probe approximation, the bulk action is renormalized through the addition of the boundary terms introduced above, giving rise to (4.62) as a result. To this we have to add the on-shell D8-brane action (g 9,M N is the induced metric on the worldvolume; its on-shell value on the U-shaped embedding is just the x 4 = const. slice of the original background geometry)  The on-shell DBI action thus has two divergent terms in the UV and needs to be renormalized.
We get the required counterterms in the probe approximation, first going to the "dual frame" metric ds 2 ∼ e −2φ/3 ds 2 [21] and then reducing on S 4 . 16 In the dual frame, the metric is asymptotically AdS, so one can figure out the counterterms as the standard volume and GH ones on that background. The resulting counterterm boundary action (already uplifted back to 8d)S can be written in terms of covariant pieces in the induced boundary metric in the dual frame. Going back to the original metric and adapting everything to the smeared case, the resulting counterterms read (4.68) 16 See also [26] for the application of the procedure for other probe branes.
Evaluating this explicitly one gets In order to cancel the divergences in the probe D8 action (4.66) we have to choose Factors in g s , l s and R appearing in the counterterms are remnants of the common asymptotic conditions (e.g. on the flux of F 4 ) for the backgrounds for which these counterterms have to be added. For the reader convenience, the counterterm for the D8-branes in standard eight dimensional form reads 17 All in all, the renormalized bulk+flavor brane action in the probe approximation reads Accordingly, the free energy density (and hence the pressure p) is (using (2.2) and (3.13)) (4.73)

The D8 counterterms in the backreacted case
In the backreacted case, the D8 counterterms introduced above are not enough to cancel the divergences, due to extra terms in the bulk gravity fields. A counterterm action of the form (4.68) is still needed, but the coefficients m 1 , m 2 given in (4.70) have to be replaced by . (4.74) One can equivalently argue that in the backreacted case one has to add to (4.69) novel counterterms with the same structure and novel coefficients. All in all the backreacted renormalized on-shell action reads (4.75) 17 Remember that in our notation T 8 = (2π) −8 α −9/2 .
In order to express this in terms of field theory quantities, one has to take into account (3.33) when rewriting aV. This way one finds (4.76) Notice that this is independent of C 2 . In order for this expression to match with that obtained in the probe approximation we need to impose (4.77) Note that condition (3.31) is consistent with this constraint. As in the unflavored case, the free energy of the system is independent of the temperature for any T in the confined phase. Accordingly, the entropy density in that phase is zero 18 as it is also evident from the absence of an event horizon in the metric. Thus, the pressure p = −f equals minus the energy density: ε = f = −p. This can be also verified by explicitly computing the energy using the standard holographic rules which equate it to the (renormalized) ADM energy of the background. The non-renormalized expression reads . (4.78) Here the metric components are written in the Einstein frame ds 2 (E) = e −φ/2 ds 2 , the integral is done over the eight-dimensional slice at constant t and r = , and K 9 is the trace of the extrinsic curvature on that slice written in terms of the nine-dimensional metric on the t = const. slice. Evaluating the above expressions on our background we get Thus, the non-renormalized energy is precisely equal to the non-renormalized free energy, consistently with our expectations.

The backreacted solution in the deconfined phase
In the unflavored case [1] the deconfined phase of the theory is holographically dual to a black hole solution. The Euclideanized temporal circle shrinks at a radial position u = u T while the x 4 circle does not shrink any more. When the D8-branes have antipodal positions on the x 4 circle, the energetically favored configuration is that of a parallel stack of D8 and anti-D8 branes, such that chiral symmetry is restored [3].

Action and ansatz
We consider smearing the two stacks of branes on the transverse circle and we turn on a U (1) gauge field on the branes, which realizes a (diagonal) U (1) B symmetry in the dual field theory. The relevant action (with the flavor brane embedding put on-shell) reads The string frame metric ansatz reads now ds 2 = −e 2λ dt 2 + e 2λ dx a dx a + e 2λs dx 2 4 + l 2 s e −2ϕ dρ 2 + l 2 s e 2ν dΩ 2 4 , We adopt an electric ansatz for the gauge field (F = dA) As usual the asymptotic value of the gauge field provides the chemical potential for the dual theory. Notice that since we look for a black hole solution where the Euclideanized temporal circle shrinks, a trivial solution of the form A t = const. is not allowed any more. Instead, we have to require that the gauge field is zero at the black hole horizon. Implementation of the ansatz above (plus the usual ansatz in (2.1) for F 4 ) gives the following 1d action which has to be supported by the zero-energy constraint (to be written with care, since the kinetic term for the gauge field is now contained in V ).

Equations of motion
The equation of motion for the gauge field can be readily solved, obtaininġ where the dimensionless constant n will be related to the charge density. The remaining equations of motion (where again we reinsert the dilaton) then read Since we work at first order in N f /N c , when plugging back the gauge field into the equations of motion for the remaining functions, it is sufficient to evaluate it on the zero-th order unflavored solution. Introducing the radial coordinate we can write the zero-th order unflavored solutions as As we have done in section 3, we solve the equations above in a perturbative expansion, defining for each field where we have used the zero-th order relations between u 0 , u T and M KK , T . 19 The above definition is such that f T = f at the phase transition. Moreover, it suggests a definition of a "running" flavor coupling f (u) ∼ √ u ∼ E where E is the field theory RG energy scale. 20 The equations of motion we need to solve are thus (derivatives are w.r.t. r) . (5.14) In this paper we do not attempt solving these equations exactly in q. Instead we focus on the small charge case, keeping only the leading q 2 terms in an expansion.

(5.17)
In the small-q case, the zero-energy condition reads This condition can be satisfied if The near-horizon behavior of the functions reads (here x = e −3r/2 → 0 near the horizon) Regularity at the tip of the (t, r) Euclideanized cigar leads to the constraints which automatically fulfill the zero-energy constraint (5.20). These conditions can be further reinforced by requiring all the logarithmically divergent terms to disappear by imposing 21 b 3 = 0 . (5.23) 21 As in the confined case, this is not strictly necessary.
The UV behavior (x → 1) of the functions is given by Requiring the leading divergences to be absent amounts on having The constant terms can be eliminated by further imposing As before we work in the more general case where the conditions (5.22), (5.25) hold. We will sometimes specialize to the more restricted cases where (5.23), (5.26) are also satisfied.

Thermodynamics
In this section we present the thermodynamic of the theory computed from the gravitational background of section 5. Let us begin by first relating our dimensionless parameter q with the quark chemical potential and density. Using the holographic relation between the asymptotic value of the original field A t and the chemical potential µ we get, from (5.15), to leading order in the N f /N c expansion where we have defined the running 't Hooft coupling The quark charge density is obtained from the holographic relation δL (5) δF tρ , (6.29) where L (5) is the five dimensional Lagrangian density (obtained from the original 10d smeared one after reduction over x 4 and S 4 ) andρ = l s ρ. An explicit computation gives Notice that, as usual, the baryon charge density is defined as From expressions (6.27,6.30) it follows that when performing calculations in the grandcanonical ensemble, where the chemical potential is held fixed, we should consider that q ∼ T −2 . Instead, in the canonical ensemble q ∼ T −5 .
Requiring the holographic computation of the thermodynamic observables (using the renormalized on-shell action for the pressure and the renormalized ADM energy for the energy density) to give consistent thermodynamic relations imposes the constraint 5a 3 = a 1 + a 2 . (6.32) In the following presentation this additional constraint is implemented. For the Euclideanized black hole metric to be regular at the horizon, the relation Note that in the formulas above M KK comes directly from the length of the S x 4 circle, so it represents the flavored energy scale.
Passing to the grand-canonical ensemble, we first notice that so that the Gibbs free energy density reads where p is the pressure. The above expression satisfies the thermodynamic relations ω = f − µn q and s = −∂ω/∂T (at fixed µ). Moreover it satisfies The trace of the stress energy tensor reads The heat capacity at fixed chemical potential is so that the squared speed of sound reads Note that the dependence on the integration constants dropped out of the thermodynamical observables, as expected. The formulas above agree with those which can be obtained in the probe approximation, see e.g. [4,27], though a covariant treatment of the holographic renormalization procedure needed in that case was missed in those works. The holographic computation of the free energy requires a careful treatment of the boundary (counter)-terms. Actually these terms are the same as in section 4, with the only difference being on the coefficients.

Holographic renormalization
As for the confined solution, the on-shell action (5.5) (and so the Gibbs free energy 22 Ω) in the deconfined case is a divergent quantity and must be renormalized by appropriate 22 The Helmholtz free energy F is holographically related to the Legendre transformed on-shell action, which, once reduced to a radial integral is defined as in (5.5) with L 1d replaced by its Legendre transform L 1d = L 1d − (δL 1d /δA t )A t , where the last term has to be evaluated on the solution of the equation of motion for the gauge field. counterterms S ren E = (S E + S GH ) + S bulk c.t. + S D8 c.t. , (6.44) where the functional form of each term is as in formulas (5.1), (4.53), (4.56), (4.68) respectively. On the deconfined solution they read where now As in the confined case, we see that the "bulk counterterm" S bulk c.t. only cancels the O( 0 f T ) divergences as in the unflavored case. In order to cancel the remaining divergences, we need to introduce additional counterterms related to the D8-branes. Another possibility, often considered in the literature, would be to cancel the divergences by subtracting to the on-shell value of S E + S GH , the value of the same combination on some reference background. In our case a natural choice would be keeping as reference background the one corresponding to the confined phase. However, by computing the difference between S E + S GH here and the related on-shell value obtained in section 4 one discovers that the O( f ) divergences do not cancel. 23 This is true also in the probe approximation: the divergences related to the on-shell D8-brane actions do not cancel among deconfined and confined solutions. This is the main reason why, in this work, we perform renormalization by adding counterterms separately in the two phases.
In the probe approximation the on-shell value of the DBI action for the D8-branes in the small-q limit reads It is relevant to notice that there are no charge-dependent divergences. This is actually true for any q (the corresponding on-shell action for the probe branes can be easily computed in this case) and it was already noticed in the past (see e.g. [27]). Thus the covariant counterterms we have introduced are enough to cancel the divergences in the charged case too. Notice instead that the divergences do not cancel when (carefully) subtracting the on-shell D8-actions (6.47), (4.65) in the deconfined and confined phases. All in all, when flavors are added in the probe approximation, the additional counterterm is the one in (4.68) with The counterterm for the D8-branes in standard eight dimensional form reads 24 In the backreacted case, the divergences are canceled provided . (6.50) Using (6.33), we get that, once the relation (6.32) is implemented, the Gibbs free energy density reads which precisely matches with (6.39). The divergences appearing in the computation of the ADM energy are exactly the ones above. Thus, the counterterms needed to renormalize it are the same.

Phase diagrams
In this section we discuss the phase diagrams of the holographic model, comparing them to the lattice results for QCD. We will first focus on the finite temperature, finite baryon chemical potential setup, and then we will also consider the model at imaginary chemical potential and at finite θ angle.

The critical temperature
We want to quantify the effects of the flavor fields on the critical temperature T c at which the first order transition between the deconfined and the confined phase happens. To get T c we just need to solve the equation p conf = p deconf (T c ) which, using (4.73) and (6.39) gives, to first order in N f /N c where we have chosen to work at fixed quark chemical potential µ, using formula (6.27). It is interesting to notice that, for any N f , the contribution of the baryon chemical potential is quadratic and it is such that T c decreases when µ increases. This is in agreement with lattice QCD results. Instead, since one might naively conclude that, at zero chemical potential, the effect of the flavors is that of increasing the critical temperature (see also [24]). Of course this is a scheme-dependent statement which means that, as in QCD [28], it depends on which appropriately chosen physical observable (e.g. the string tension or the ρ-meson mass) or suitable UV mass scale (as in [29] or [30]) is held fixed when comparing theories with different numbers of flavors. In fact, while the chemical potential is a parameter of a given theory and so one can unambiguously compare the theory at different values of µ, this is not the case for N f . By varying the number of flavors one is changing the theory: comparison of two theories requires a scheme and there is no a-priori "correct scheme" -it depends on the different physical effects one wants to investigate. See [22,31,15] for some examples and discussions in the holographic context. In the present case the behavior of T c with N f changes qualitatively by changing comparison scheme. For example, one of the comparison schemes used in lattice literature consists in fixing the coupling at a certain energy scale. In our case, this produces a shift in u 0 , giving an additional term in (7.52) which depends on the solution's parameters -the sources of gauge invariant operators. Now, for some choices of these parameters, the overall effect of the flavors is still that of increasing the critical temperature. But for large enough values of parameters the effect is the opposite, provided the energy at which we fix the coupling is appropriately chosen.
Analogously, while the comparison scheme which fixes the mass of the ρ meson gives an increasing T c with N f , one can find a simple scheme where T c decreases with N f , by fixing the ratio of the string tension with M 2 KK . Thus, not having a particular reason to prefer a comparison scheme w.r.t. another one in the present discussion, we will not consider the dependence of T c on N f furthermore in the following discussion of the phase diagrams.

Finite θ and imaginary chemical potential
An interesting parallel between the phase diagram of (quenched) QCD at finite θ B ≡ µ B /T c (0) (where µ B is the imaginary baryon chemical potential 25 and T c (0) the critical temperature at µ B = 0) and that of the pure Yang-Mills theory at finite θ angle has been traced in [32]. Let us describe this relation and its realization in the holographic model at hand.
When θ is turned on, the lattice QCD partition function displays imaginary terms which prevent Monte Carlo simulations to numerically converge (see e.g. [33] for a review). This "sign problem" (analogous to the one encountered at finite, real, baryon chemical potential) can be avoided if the theory is formulated at imaginary values of θ. In recent years lattice studies of the QCD phase diagram at finite imaginary chemical potential [34] or θ parameter [35] have been performed. Assuming the (free energy of the) theory to be analytic around θ B = 0 or θ = 0 one can perform an analytic continuation towards real values (and a further continuation to the continuum limit) to get information on the phase diagram of the theory at small real values of these parameters.
Using this trick, for example, the effect of a finite small θ parameter on the deconfining temperature T c of pure SU (3) Yang-Mills was considered in [6]. The result was found to be in qualitative agreement with the following formula give , with R θ ∼ 0.0175, while large N c arguments suggest r θ ∼ 0.25 [6]. In both cases, the deconfining temperature is found to be a quadratically decreasing function of θ.
Going beyond the small θ regime is not possible on the lattice at the moment, due to the above mentioned sign problem. However, at least in the large-N c limit of Yang-Mills, it is possible to figure out what is the structure of θ-dependent observables just by using simple scaling arguments. In 't Hooft's large N c limit, θ-dependence can remain non trivial provided θ/N c is held fixed. In order to reconcile this condition with, say, the periodicity condition ε(θ) = ε(θ + 2πk) on the vacuum energy density of the theory, Witten proposed that the latter (which in turn should have an absolute minimum at θ = 0) should have the form [36,37] where b is a model-dependent dimensional coefficient (independent on N c ) and k is an integer. The vacuum energy is thus expected to have discontinuities at θ = (2k + 1)π. At each of these points (where CP invariance is recovered), a first order quantum phase transition occurs between two different phases of the theory where CP-symmetry is broken. While the above structure cannot be verified yet on the lattice, formula (7.55) found an explicit realization in [37]. There, the computation was done, using the holographic gauge/gravity duality, for the unflavored version [1] of the model we have focused on in this work. As it has also been argued in recent literature (see e.g. [32]), at finite temperature the same structure as in (7.55) should be inherited by the Yang-Mills free energy density, provided the theory stays in the confined phase. The (T, θ) phase diagram should thus present vertical first order lines departing from the (T = 0, θ = (2k + 1)π) points. These lines are however expected to end at around T = T c (0), where a parabolic phase boundary described by an equation like (7.54) should be met. The resulting phase diagram will thus present itself as an "arcade", with first order "columns" emerging from (T = 0, θ = (2k + 1)π) ending at tricritical points at T c = T c (θ = (2k + 1)π) where they bifurcate in two first order "arcs". See figure 2 for a sketch.
As observed in [6,32], this structure resembles that of an inverted Roberge-Weiss phase diagram [38]. The latter describes (quenched) QCD at finite imaginary baryon chemical potential µ B . In this case the high temperature free energy (which can be accessed perturbatively) shows discontinuities at θ B = (2k + 1)π. For an SU (N c ) theory with N f massless Weyl fermions in the fundamental representation, the tree level result for the θ B -dependent part of the free energy density f (θ B ) is The phase diagram (sketched in figure 2) is then expected to exhibit first order "columns" emerging from the points (T → ∞, θ B = (2k + 1)π). These lines should end on the phase boundary at T ∼ T c (θ B ) which is shown to have a quadratic dependence on θ B on the lattice (with T c (θ B ) increasing with θ 2 B at small θ B ). The whole phase diagram looks like an inverted arcade with arcs around T c (0).
In the following we will show that the holographic QCD model studied in this paper gives a concrete realization of the structure of the above described (T, θ B ) and (T, θ) phase diagrams, which are, then, qualitatively "dually" related as proposed in [32].

Holographic QCD at finite T and θ B
In the case of imaginary baryon chemical potential the result (7.52) readily gives Thus, the critical temperature quadratically increases with θ B . This behavior, for the Witten-Sakai-Sugimoto model, has been already found in [39]. In the latter paper, generalizing the results of [40], it has also been shown that the free energy has first order discontinuities at θ B = (2k + 1)π, so that the phase diagram is precisely as in figure 2.

Holographic Yang-Mills at finite T and θ
We now want to construct the phase diagram of the holographic model at finite temperature and finite θ-angle. This requires a few considerations, whose outcome is formula (7.64) -the reader uninterested in the details can skip the derivation. 26 The discussion concerns the unfavored theory dual to the backgrounds of section 2. Let us take 0 ≤ θ < π and start by reviewing the known results at T = 0 [37], thus considering the background (2.1). The θ-angle is identified with the integral of a Ramond-Ramond one-form where "cig" is the (u, x 4 ) cigar, F 2 = dC 1 = F u,x 4 du∧dx 4 and we have used Stokes's theorem.
Since the S x 4 circle shrinks to zero size at u = u 0 , regularity imposes C(u 0 ) = 0. With this condition, the solution of the equation of motion for C 1 gives [37,43] Remember that β 4 is the length of the S x 4 circle. Under the holographic gauge/gravity map, the θ-dependent contribution to the field theory vacuum energy density is mapped into the part of the gravity action (in Euclidean signature) involving F 2 S[F 2 ] = − 1 4π(2πl s ) 6 dx 10 √ −g|F 2 | 2 , (7.60) 26 See also [42] for relevant related results.
which has to be computed on-shell on the solution for F 2 found above. Using the holographic relation e −V 4 ε(θ) ≈ e −S[F 2 ] , where V 4 is the (infinite) 4d Euclidean volume, one gets (see [2]) This formula is thus of the anticipated form given in (7.55), restricted to the 0 ≤ θ < π range.
The T < T c geometry is precisely the same as the T = 0 one, but for the compactification of the Euclidean time direction. Therefore, the holographic relation e −βF ≈ e −S between the free energy F and the on-shell Euclidean gravity action, implies that the free energy density f (θ) in the confined phase is given by f (θ) = ε(θ), with ε(θ) given in (7.61). In particular this means that f (θ), in this phase, is independent of T and it has discontinuities at θ = (2k + 1)π exactly as in the T = 0 case. The phase diagram at T < T c will thus show first order phase transition vertical lines emerging from the (T = 0, θ = (2k + 1)π) points.
At T > T c , instead, the free energy of the model does not depend on θ (this has to be read as the leading effect in the large N c limit), since a trivial solution for F 2 is now admitted. In fact (see the background (2.3)), the (u, x 4 ) subspace is a now cylinder, i.e. the S x 4 circle does not shrink anymore. Such a geometry allows for a trivial solution C(u) ∼ θ, i.e. F 2 = 0 of the equation of motion of C 1 . This is actually the solution which minimizes the gravity action. The latter (which depends on F 2 ) will thus be independent of θ. As a result, the field theory free energy will not depend on the θ angle. We now show that this implies a structure of the form given in (7.54) for the θ-dependence of the deconfining temperature.
When θ = 0, the difference between the free energy density at T < T c and that at T > T c in the unflavored case is given by the difference of (6.39) and (4.73) with f = f T = 0 (see [3] 27 ). If we want to compute the same quantity at finite θ it just suffices to add the T < T c contribution f (θ) = ε(θ) as given in (7.61). The result is (7.63) Thus, the critical temperature for the phase transition at fixed θ is found to be given by where T c (0) = M KK /2π is the value at θ = 0, and we have expanded the parenthesis considering that the second term is suppressed w.r.t. 1 for small θ/N c .
Thus, T c is a quadratically decreasing function of θ as in (7.54) with coefficient The whole result is thus in qualitative agreement with large-N c considerations and lattice Yang-Mills results obtained by analytic continuation from imaginary θ [6]. Furthermore, the periodicity of the physics with θ dictates the whole phase diagram to have the periodic structure of an inverted Roberge-Weiss [38] one (as in figure 2), with cusps joining the first-order lines emanating from T = 0 at θ = (2n + 1)π. All in all, the structure of the phase diagram of the holographic model is the same one advocated in lattice QCD ( fig. 2), both at finite θ and imaginary chemical potential. The qualitative "duality" 28 of the two phase diagrams in the holographic model depends crucially on the symmetry of the T < T c and T > T c backgrounds (2.1), (2.3) under the exchange of the S x 0 and S x 4 circles. In fact, the latter determine the behavior of the free energy density with θ and θ B , as they support the gravity one-forms C 1 and A 1 dual to these quantities.
The quadratic dependence of f on θ and θ B just follows from the Maxwell form of the action for the one-form fields. 29 As usual, the holographic description nicely geometrizes properties of the dual field theory.

Degrees of freedom for the phase diagrams
In the stringy picture there are degrees of freedom which emerge naturally in the confined and deconfined phases of the theory. In fact, D0-branes are electrically charged under C 1 (related to the θ angle). This provides a way for identifying the instantons in the model as Euclidean D0-branes wrapped along the S x 4 circle (see figure 3). 30 Since at T > T c the circle does not shrink, this configuration is stable and maps to a instanton gas configuration in the dual field theory. At T < T c this picture does not hold anymore and this is realized by the fact that the wrapped D0-brane is now shrinking with the cycle.
Analogously, we propose that the degrees of freedom replacing the instantons in the confined phase are baryon vertices winding around the Euclidean S x 0 time circle. This idea is supported by the fact that a baryon vertex, which is identified with aD0-brane obtained by wrapping D4-branes on the transverse S 4 , is electrically charged under the U (1) B gauge field A t . In the confined phase the fact that the (x 0 , u) subspace is a cylinder, allows for the wrappedD0-brane configuration to be stable (see figure 3). This suggests a possible free "baryon gas" interpretation in the dual field theory. When T > T c this configuration is instead unstable and the above interpretation cannot be used anymore. 28 Needless to say, this is not a precise duality but only an interesting qualitative analogy. 29 This is true for A 1 in the small charge regime we are considering. At larger charges the non linear terms in A 2 1 in the Born-Infeld action could correct this behavior. 30 See e.g. [41,42] and references therein. Therefore, the D0-brane in the (T, θ) case and theD0-brane in the (T, θ B ) case, play a "dual" role. Again, this is due to the S x 0 ↔ S x 4 symmetry of the backgrounds, together with the existence of two types of "effective D0-branes" in the theory.

Conclusions and discussion
In this paper we have constructed gravity backgrounds corresponding to the Witten-Sakai-Sugimoto (WSS) model with backreaction of the D8-branes, that is with dual dynamical flavors. We have provided solutions at small (or zero) and large temperature and at finite chemical potential (or charge density). The solutions are suitable for the study of the influence of dynamical flavors on physical observables in a setting (the WSS theory) which represents the top-down model closest to (planar) QCD.
In particular, in this paper we have concentrated our attention on the phase diagrams, re-deriving the thermodynamic of the system and discussing it with an eye to recent lattice QCD literature. We have observed qualitative similarities between the QCD and WSS theories in the finite temperature, finite baryon chemical potential regimes. The similarities extend to the finite θ angle and finite imaginary chemical potential regimes too. As usual, the holographic approach, while describing a theory which is not exactly (planar) QCD, has the virtue of geometrizing certain aspects of the physics and suggesting information on the possible relevant degrees of freedom.
Moreover, we have provided covariant counterterms for the D8-branes in the Witten-Sakai-Sugimoto model allowing for its standard holographic renormalization. As far as we know, these counterterms have never been explicitly reported in the vast literature on the model. Our main motivation for deriving these counterterms has been to holographically renormalize the free energies, in order to study the phase diagrams and the confinement/deconfinement transition. In many examples of holographic theories (see e.g. [44,31] for flavored cases) holographic renormalization is not needed, since the background subtraction method gives a finite result for the difference of free energies. Interestingly, this is not the case for the WSS model. In fact, the D8-branes have different sub-leading divergences in the T < T c and T > T c phases, so they need to be renormalized independently (this is the reason for having different coefficients of the counterterms in the two cases). Note also that these subleading divergences carry information on the IR region of the theory (they depend on u 0 or u T ), suggesting that there could be no clear-cut separation of the UV and IR regimes.
Anyway, despite the presence of a Landau pole due to irrelevant operators, the theory admits a standard holographic renormalization. This might depend on the fact that we have worked at leading order in the flavor counting parameter f ∼ λ 4 N f /N c , as advocated in [45]. The situation might be more involved at higher orders. Moreover, while we have renormalized the gravity action, relevant for the computation of correlators of local operators, a separate renormalization prescription is needed for non-local observables, such as Wilson loops. We plan to come back to this interesting issue in the future.
As said, the solutions we have constructed in this paper include the backreaction of the D8-branes at leading order in the flavor counting parameter f . This allows to find analytic solutions. Going beyond the leading order would face two type of problems. The first one is technical: finding analytic solution at higher orders appears to be a hopeless task and one should resort to a numeric investigation, which is not guaranteed to be straightforward. In fact, higher order terms in f require to add gravity fields in the setting, especially in the charged case (as in the D3 − D7 system [31,15,46,47,16,17]). This causes the proliferation of degrees of freedom and a complication of the equations to be solved. Already finding analytic IR asymptotics can be problematic in this situation. Moreover, the presence of the UV Landau pole demands to solve the equations up to a finite cut-off, which brings in the game a certain degree of uncertainty and technical complications. Finally, at higher order in f , there could be corrections to the DBI action, such as blackfold terms, taking into account thermalization of the brane degrees of freedom [48]; it would be extremely difficult to include such contributions in the WSS theory.
On top of this technical complications, there is a second type of problem, having to do with the energy scale of the Landau pole. As one can readily realize, and as can be seen from (3.36), the energy scale of the Landau pole is inversely proportional to f . Requiring a finite separation between the Landau pole scale and the IR physics, forces us to stay in the small f regime. Thus, the leading f -solutions we have constructed furnish the main perturbative contribution to the exact solution. 31 In this sense, knowing the exact solution is not a urgent priority.
Since the small f regime is mandatory in this model, one could think that the probe approximation is sufficient to capture all the relevant physics. This is only partially true: the probe and backreaction approaches are more complementary than equivalent. In fact, whenever calculations can be performed both in the probe approximation and in the leading f case, the results should coincide. Thus, one would likely prefer to use the easier approach of the probe approximation. On the other hand, for many observables one does not have a given prescription on how to extract information from the probe -this has to be derived case by case. Instead, with the backreacted solution one can use the standard prescriptions derived for the unflavored theory, involving only the background, to compute observables in a straightforward way (we have given examples in section 3.1).
Coming to our other approximation, it would be certainly useful and worthwhile to try and solve equations (5.14) beyond the small charge regime. There might exist analytic solutions, which would allow to study the physics in a regime where lattice QCD simulations are problematic. Other obvious venues for future research can include, among the others, the study of transport properties of the theory (as done in [49] for the D3 − D7 system) and the flavor and charge dependence of the entanglement entropy [50].
Moreover, it would be interesting to consider the backreaction of explicit baryonic sources, so to extend the results of, e.g. [5,51].