Non-equilibrium physics at a holographic chiral phase transition

The D3/D7 system holographically describes an N=2 gauge theory which spontaneously breaks a chiral symmetry by the formation of a quark condensate in the presence of a magnetic field. At finite temperature it displays a first order phase transition. We study out of equilibrium dynamics associated with this transition by placing probe D7 branes in a geometry describing a boost-invariant expanding or contracting plasma. We use an adiabatic approximation to track the evolution of the quark condensate in a heated system and reproduce the phase structure expected from equilibrium dynamics. We then study solutions of the full partial differential equation that describes the evolution of out of equilibrium configurations to provide a complete description of the phase transition including describing aspects of bubble formation.


Introduction
Thermal phase transitions are a crucial aspect of the evolution of the Universe after the Big Bang and also in the physics of heavy ion collisions. We have traditionally lacked tools to study these transitions though in strongly coupled systems such as QCD. The AdS/CFT correspondence [1][2][3], which gives a weakly coupled string/gravity description of a class of strongly coupled gauge theories, offers the chance to study similar transitions in detail.
In this paper we study a simple dual of a theory with gauge fields and quarks which has (in the presence of a magnetic field) a first order chiral symmetry restoring phase transition as it is heated. Previous analysis of this gauge theory (using probe branes in an AdS-Schwarzschild black hole geometry) has explored the first order transition for equilibrium (time independent) configurations. The heating or cooling of the system can be studied thanks to the boost-invariant expanding or contracting plasma geometry of Janik and Peschanski [4]. That geometry, which we will review below, has a moving black hole horizon describing the changing temperature in the gauge theory. We will study probe branes in this geometry to learn more about the first order phase transition out of equilibrium.
The particular duality we concentrate on is the simplest example of holography with fundamental quark fields [5][6][7][8][9][10]. We do not consider the specific degrees of freedom of the theory too crucial -it is some strongly coupled gauge theory that displays a chiral phase JHEP01(2011)050 transition. We hope, in the spirit of AdS/QCD models [11,12], that it reflects broad aspects of many strongly coupled systems. The specific gauge theory is constructed from the D3/D7 system in type IIB string theory which we will describe further below. The theory is the large N c N = 4 U(N c ) gauge theory with a small number of quark hypermultiplets. We will work in the quenched approximation [7] (appropriate when N f ≪ N c ) which on the gravity dual side corresponds to treating the D7 branes as probes in the metric generated by the D3 branes. There is a U(1) chiral symmetry (a remnant of the SU(4) R-symmetry of the N = 4 theory) which is broken when a quark condensate forms [9,13,14]. Several mechanisms for triggering this condensation have been explored. The cleanest is when a background magnetic field is introduced [15][16][17][18][19][20][21][22][23] so we will use that mechanism here. The reader might want to loosely view the B field as simply the introduction of a conformal symmetry and supersymmetry breaking parameter that triggers the strong dynamics to cause the symmetry breaking. Physically, magnetic fields may be strong during structure formation in the early universe, in particular during the epoch of the QCD primordial phase transition, and in non-central high-energy heavy ion collisions studied at RHIC. Running of the coupling in the holographic dual also causes quark condensation as has been shown in back-reacted dilaton flow geometries [9,24] and models with a phenomenologically imposed dilaton profile [25]. The quark condensate can be determined in these models and an effective IR quark mass is generated. The theories display a massless pion-like Goldstone field and a massive sigma field (since we are at large N c it is stable) that is the effective Higgs particle.
The equilibrium finite temperature behaviour of the theory with a magnetic field has been studied in [15-17, 22, 23]. Finite temperature can be included through the presence of an AdS Schwarzschild black hole. At a critical temperature the D7 embedding flips from a chiral symmetry breaking embedding away from the horizon to a symmetry preserving embedding that enters the horizon. The transition is therefore also associated with meson melting [26][27][28] -for embeddings away from the horizon there are regular linearized fluctuations describing the meson spectrum of the theory. For the embedding that enters the horizon there are only in-falling quasi-normal modes describing unstable plasma fluctuations. In terms of the quarks of the theory the high temperature phase is analogous to the quark gluon plasma phase in QCD whilst the low energy phase is more akin to the hadronic phase of QCD. It should be noted though that the gluonic degrees of freedom deconfine at any finite temperature so the analogy is imperfect.
The crucial extra ingredient we shall add to this story in this paper is provided by the boost-invariant expanding or contracting plasma geometry of Janik and Peschanski [4]. This geometry has a black hole whose horizon moves away from the boundary in time as the N = 4 plasma it describes expands and cools. The time reversed solution at zero viscosity describes a heating plasma and we will find it useful to discuss that scenario too below. The geometry is a late time expansion (when the black hole is small) in powers of inverse time. However, by controlling the strength of the magnetic field on the D7 probes, felt by the quarks, we can arrange to place the chiral phase transition at any point in the evolution so the expansion is sufficient to fully study the transition.

JHEP01(2011)050
We will place a D7 brane 1 into the expanding plasma geometry and determine the partial differential equation (PDE) that describes the time dependence of its embedding. A good first approximation to the transition behaviour is provided by the equilibrium results with the temperature replaced by the temperature as a function of time from the moving background. In fact we will see that that is an extremely good approximation when talking about the slow or adiabatic heating or cooling of vacuum configurations at temperatures even close to the phase transition. If the chiral symmetry breaking embedding is heated though, the local minima in the effective potential associated with that embedding is eventually lost and the configuration becomes an out of equilibrium configuration. The equilibrium results can not describe the subsequent evolution. Similarly excited vacuum states must be studied through the full PDE system.
We first turn to an approximation to the PDE. The solution can be power expanded in inverse powers of time. This reduces the PDE to a system of ordinary differential equations (ODEs) that are much easier to solve. This analysis was first done in [31] where the time evolution of the high temperature phase of the pure N = 2 theory was studied. Using this technique we solve the ODE system and find the moving D7 solution. This method assumes that the initial vacuum state is exactly a maximum or minimum of the zero temperature effective potential. Again the heating is essentially adiabatic in nature. The expansion breaks down if the super-heated state ceases to be an extrema of the effective potential. This method allows us to confirm the success of the equilibrium derived results although in fact the expansion breaks down before the equilibrium results deviate from the full PDE solutions.
The most interesting out of equilibrium questions lie beyond the adiabatic approximation though. In a physical first order transition quantities such as the condensate do not jump but the vacuum state instead performs a fast roll from one configuration to another. The timing of that transition can and most likely will be spatially dependent i.e. bubbles of the true vacuum will form and grow. We turn to solving the full PDEs to study these phenomena (with care this can be done using in built PDE solvers in for example Mathematica). For example we are able to watch the super-heated chiral symmetry breaking phase roll to the symmetric phase. We can also simulate initial conditions with some extra energy (which might for example come from thermal fluctuations) and see such configurations transition away from the super-heated vacuum before the meta-stable vacuum has disappeared in the effective potential. This allows us to confirm some elements of the transition such as the length of time in which there is a mixed phase. We also are effectively watching a very large bubble form.
Our main result then is to have developed numerical techniques that let us reproduce the phase structure using the PDE solutions and to describe non-equilibrium physics that is necessarily present in the first order transition.
In our final section we also analyze the ODE expansion approximation to the PDE solutions for the D7 embedding to make clear the full dependence of the solutions on the JHEP01(2011)050 magnetic field value. In that case the dependence is available analytically. We also show the effects of the viscosity.

Holographic descriptions
In this section we review the N = 2 gauge theory with a magnetic field and its holographic description. We discuss the theory's finite temperature chiral phase transition in an equilibrium description. We then review how to study flavour physics in a nonequilibrium set-up using holography.

Chiral symmetry transition
We will study the D3/D7 brane model at finite temperature and with a magnetic field [15-17, 22, 23]. The magnetic field, which causes chiral symmetry breaking, competes with the temperature that prefers to restore chiral symmetry. Consequently there is a first order phase transition.
The N = 4 SU(N ) gauge theory at finite temperature has a holographic description in terms of an AdS 5 black hole geometry which can be written as where φ is a U(1) angle and w = ρ 2 + L 2 , ρ := w sin φ, L := w cos φ and , Note R 4 = 4πg s N α ′2 , and the temperature is given by w H := πR 2 T .
Quenched (N f ≪ N ) N =2 quark superfields can be included in the N = 4 SU(N ) gauge theory through probe D7 branes in the geometry [7]. The D7 probe can be described by its Dirac-Born-Infeld (DBI) action where P [G] ab is the pullback of the metric and F ab is the gauge field living on the D7 world volume. We will use F ab to introduce a constant magnetic field [15][16][17], We embed the D7 brane assuming only ρ dependence: L(ρ) at constant φ. The full DBI action we will consider is then where ǫ 3 is a volume element on the 3-sphere and  with the dimensionless variables defined as Note that we can use the magnetic field value as the intrinsic scale of conformal symmetry breaking in the theory -that is we can rescale L and ρ by B. The Euclidean on-shell Lagrangian (−L) is interpreted as the free energy density ( Ω). In all cases the embeddings become flat at large ρ taking the form where m and c are identified with the quark mass and the quark condensate, respectively. Since we are interested in spontaneous symmetry breaking we impose m = 0 in the UV ( ρ → ∞). Then the value of c is determined by requiring regularity in the IR ( ρ ≪ 1). If there is more than one solution then we choose the one giving the minimum free energy. The results are displayed in figure 1. At low temperatures T ≪ B, the black hole is small and the embeddings are repelled from the origin of the L − ρ plane (figure 1 (left)). This behaviour is a result of the inverse powers of w 4 , when w H ≪ 1, in the last term of the Lagrangian (2.6), which causes the action to grow if the D7 approaches the origin. Consequently c is non-zero and chiral symmetry is broken. If the temperature is allowed to rise sufficiently then the black hole horizon grows to mask the area of the plane in which the inverse w 4 terms in the Lagrangian are large. At a critical value of T

JHEP01(2011)050
(w H = 0.2516, figure 1 (middle)) the benefit to the m = 0 embedding of curving off the axis becomes disfavoured and it instead lies along theρ axis -chiral symmetry breaking switches off. This transition is first order since the condensate vanishes discontinuously, which corresponds to the embedding change from the blue to the red one. At higher temperatures (figure 1 (right)), the embedding stays flat, c = 0, as expected in the chirally symmetric phase. The corresponding condensate vs. temperature (−c-w H ) plot is shown in figure 1 (bottom).

A boost-invariant expanding plasma
The geometry (2.1) is dual to a system in thermodynamical equilibrium and therefore not suitable for the description of the chiral phase transition in a rapidly expanding plasma, in which the transition is basically a non-equilibrium process. Boost-invariant expanding N = 4 SU(N ) plasmas out of equilibrium can however be described by the time-dependent background found in [4] (see [32,33] for a review). In the following we review the basic features of this background and discuss the embedding of probe D7 branes dual to quenched flavours ("quarks") in the plasma.
The boost-invariant geometry is a 5D spacetime with coordinates {τ, y, x ⊥ (= x 1 , x 2 ), z}, which, apart from the holographic direction z, parameterise the 4D spacetime on the boundary. The longitudinal position plane is parameterised by the proper time τ and rapidity y (related to x 0,3 as x 0 = τ cosh y, x 3 = τ sinh y), the transverse coordinates are collected in x ⊥ . We also add a five-sphere to obtain a full type IIB supergravity background. The metric is then of the form where R is the radius of the AdS 5 space. At late times, the coefficients can be expanded to first order as [4,34,35] a(τ, z) = ln ε 0 is a free parameter of mass dimension 8/3 and is related to the energy density, while η 0 is related to the shear viscosity. v is a scaling parameter valid at large τ . Note that a(τ, z), b(τ, z) and c(τ, z) are expanded around τ = ∞ in powers of 1/τ 2/3 as [4,34,35] a(τ, z) = ∞ n=0 a n (v) 1

JHEP01(2011)050
and similarly for b(τ, z) and c(τ, z). The coefficients are functions of a scaling parameter v only. Due to this scaling behaviour the complicated PDE Einstein equations can be reduced to an ODE system, which allows analytic solutions. At early time τ ≪ 1 there is no scaling behaviour and we should solve the full PDEs.
To find which gauge theory state corresponds to this bulk metric according to the AdS/CFT dictionary let us expand the metric around z = 0, The leading terms (of order O(z 0 )) of the metric elements (2.13) is simply the Minkowski metric in the τ −y coordinate system, and the subleading terms (of order O(z 4 )) correspond to the expectation value of the energy-momentum tensor, i.e.
This energy-momentum tensor is precisely that of a longitudinal viscous boost-invariant N = 4 SYM conformal plasma with finite η.
In g τ τ we recover a time-dependent emblackening factor, which describes a moving horizon. The size of the horizon determines the time-dependence of the 'temperature' as [31,34] (2.17) If we assume that the time-dependent entropy density s(τ ) has the same form as in the static case, s(τ ) = (π 2 /2)N 2 c T (τ ) 3 , then the ratio η/s can be computed as [35] which at large τ agrees with the known static bound. Note that the numerical value of η 0 in (2.11) is crucial to get 1/4π.
has been factored out to apply the AdS/CFT dictionary (2.14)

JHEP01(2011)050
Finally we note the uniqueness of the gravity solution. There is a potential singularity is only regular at each order if one makes the specific choices above. We must choose −1/3 for the power of τ in (2.11) to make R 0 regular [4]. We must also choose (2.10) with the specific choice of the numerical value η 0 as in (2.11) to make R 2 regular [35]. (R 1 is always regular and does not give any constraint.) Flavours in this background were first studied in [31]. Knowing the explicit background geometry (2.9) we can study flavour physics using the D7 brane DBI action. The action schematically reads where A and B are complicated but known functions of τ and ρ, see (3.5) and (3.6) in the next section. L = L(τ, ρ) is the embedding profile of the D7 brane and is assumed to be a function of τ and ρ.
The equation of motion coming from (2.19) is a non-linear PDE. However it can be semi-analytically solved by the late-time expansion This expansion inherently assumes that the late time configuration is precisely the equilibrium vacuum and when we heat it there is no excess energy. This leads to an adiabatic approximation. The fraction 1/3 is chosen because all exponents of τ in A and B are integer multiples of 1/3 and the constant m reflects the fact that the embeddings get flat at large τ . Here we only consider a large bare quark mass (m ≫ 1) and Minkowski embeddings, i.e. D7 configurations which end well above the black hole horizon. This reduces the partial differential equation to an infinite set of ordinary differential equations given by

21)
where we do not consider terms with i > 13 since they are beyond the validity regime of our approximation of the boost-invariant metric. The asymptotic solution at large ρ is We set m i = 0 so the bare quark mass is zero and not time-dependent. We also impose the condition f ′ i (0) = 0 for regularity. With these boundary conditions, we find f i = 0 except for f 8 , f 10 . As a result, one gets The condensate approaches zero as ∼ − τ −8/3 and the viscosity has a "dragging" effect ∼ + τ −10/3 . We close this section with some technical remarks on black hole embeddings. In the Fefferman-Graham (FG) coordinate system (2.9) one cannot approach the horizon for fixed τ as in the Schwarzschild black hole metric. There is also an extended background written in terms of Eddington-Finkelstein (EF) coordinates, where the spacetime is well defined across the horizon [36,37] and black hole embeddings may be described more consistently there. However we have found that the embedding configuration is not easy to handle in those coordinates. Below, to enable us to study black hole embeddings, we use FG coordinates but with a cut-off slightly above the horizon. Since the FG coordinate system is a good patch near the horizon at large τ , by restricting ourselves to large τ , we may capture the essential physics of the embedding. This will allow us to go beyond the results in [31].

Out of equilibrium description of the chiral phase transition
In the previous section we reviewed the holographic description of dynamical flavours in an expanding plasma [31]. A boost-invariant background with embedded D7 branes is however not sufficient to study the chiral transition. In that case the D7 embeddings are always flat in the chiral limit (m → 0) and the system is always in the chiral symmetric phase with vanishing quark condensate. In order to describe the transition to the chiral broken phase, we need a repulsive effect to compete against the attractive force of the black hole. As in the static case reviewed in section 2.1, this can be achieved by turning on a magnetic field. In this section, we will therefore consider D7 branes with a world-volume magnetic field in the dual geometry of an expanding plasma. In this way we will find the holographic dual of a chiral transition and deduce the dynamic effective potential for the time-dependent quark condensate.

D7 flavour brane action
The background metric for a boost-invariant expanding plasma can be written as where r 2 ≡ ρ 2 + L 2 . The S 5 part is written as in (2.1) and the AdS 5 part follows from (2.9) with a change z → R 2 /r, i.e. a(τ, r) ≡ a(τ, z → R 2 /r) in (2.10) and similarly for b(τ, r) and c(τ, r).
We are interested in time-dependent D7 brane embeddings of the type L = L(τ, ρ). The corresponding DBI action is where we turn on a constant magnetic field [15][16][17], in order to induce the chiral symmetry breaking. More explicitly, the D7-brane action reads where A = e a/2+b/2+c , B = e −a , C = e −2c , η 0 = 1/(2 1/2 3 3/4 ) as in (2.11), and ǫ 3 is the volume form of the three-sphere. For vanishing B, the action reduces to that in (2.19).
It turns out to be convenient to work with the rescaled variables For R = 1 the action then reads with A, B, C as in (3.5)-(3.7) but now expressed in terms of τ , ε 0 , ρ, and L.

Naive equilibrium based approximation
An obvious first approximation to understanding the time dependent chiral phase transition in this set up is to use the equilibrium results from section 2.1. There we described how the D7 embedding behaved in the background of a fixed size black hole and determined the quark condensate as a function of temperature. In the cooling plasma geometry of section 2.2 the black hole horizon moves as a function of time τ as If this heating were very slow (as it is at large τ ) we would expect to be able to plot the quark condensate against τ by simply substituting for the temperature T in the equilibrium We can recast the c vs. τ plot schematically as a time dependent effective potential by fitting it at each value of τ to a potential of the form using the values of c at the extrema to fix the parameters (the overall scale is not set but the figure is intended to only be schematic). We plot this in figure 2b. This is the standard picture of a first order transition. It is important to interpret these plots correctly. Firstly note that there is a c = 0 solution for all τ . If we begin at high temperature with the symmetric flat embedding and cool slowly or adiabatically (i.e. move to large τ ) then we can remain in that embedding for all τ . Above the time τ a = 20.6 two extra solutions for the condensate develop and these, as we will discuss, trigger the first order transition. At large τ the c = 0 flat embedding becomes a local maximum of the effective potential. To stay in the flat embedding for all τ is the extreme limit of super-cooling the high temperature phase into an unstable vacuum state.
The top trajectory in the c vs. τ plot is best thought about in the time reversed solution that is heating up. At early times (large τ in the plot) the solution is the T = 0 symmetry breaking D7 embedding vacuum. Now as we heat adiabatically the condensate tracks along the trajectory to smaller τ . The solution ceases to exist at τ a , where the cusp in the (black) curve is, indicating that the minimum of the effective potential corresponding to this solution has ceased to exist at τ a . At this point the brane will move quickly, as an out of equilibrium configuration, ending as the flat embedding or oscillating about it. The real-time evolution of the configuration is represented by the path of red dots in figure 2b. Starting from the global minimum with c = 0, the minimum is lifted up and JHEP01(2011)050 the configuration moves adiabatically. At τ a the minimum disappears and the system rolls down the potential to the true vacuum with c = 0. The equilibrium results can not tell us about this motion. Now we also see the role of the middle arc in the plot -these solutions are symmetry breaking embeddings that end on the horizon and they are a local maximum of the potential lying between the two minima (the embeddings of the top trajectory and the flat embedding).
In fact the first order transition point between the symmetry preserving and the symmetry breaking embeddings can be computed from the free energy in the equilibrium computation. This transition occurs at τ c = 27.5, where the local and the global minimum in the effective potential interchange their roles, i.e. for τ c ≥ 27.5 the symmetry breaking vacuum (with c = 0) becomes the global minimum. Another important time determined by the equilibrium computation is τ b = 34.9 where the mixed phase ends at large τ .
Our first task in solving the PDEs for the D7 brane motion that result from (3.4) is to demonstrate that this description is essentially correct -we will see that it is. We will try to find confirmation of the times τ a−c of these events. The more interesting task is that we will be able to follow the evolution of a particular initial condition through the phase transition. Of course in the first order transition the brane configuration does not discontinuously leap between the symmetry breaking embedding to the flat embedding but evolves continuously. We will provide solutions for this evolution.
Another interesting phenomena associated with a first order transition is bubble formation. In real systems thermal energy will lead to volumes of space which have more energy than an equilibrium like state in a local minimum of the potential at any particular time. These volumes may "climb" over the potential hill to the other local minimum during the mixed phase period shown in figure 2. These bubbles then grow or contract triggering the phase transition around τ c ending any super-heated or cooled phase. We will not look at (x 3 ) spatial dependent brane embeddings. However, we can inject kinetic energy in the holographic directions of our description into the brane configuration before we heat or cool it to see the configurations moving more quickly between the two local minima than the lowest energy configuration. This will allow us to test the time period in which the mixed phase exists in the out of equilibrium problem.

Adiabatic dynamic D7 brane embeddings
Our first study of the PDEs describing the D7 embedding in the expanding plasma geometry will be to study adiabatic expansion. The non-linear partial differential equations resulting from (3.10) can again be transformed into a set of second-order ordinary differential equations as we reviewed for the case with no symmetry breaking in section 2.2. For that, we use the late-time expansion Thus we will be following the evolution of the symmetry breaking embedding as τ decreases.
There is also the solution f i = 0 which corresponds to the symmetric embedding. Note that, contrary to (2.20), a nontrivial asymptotic embedding f 0 ( ρ) is assumed because of the repelling effect of the magnetic field. The equations for every f i ( ρ) can be obtained order by order in τ − 1 3 , which will be solved recursively. Unlike (2.21) these equations are quite lengthy and will not be presented here. They do not allow for an analytic solution such as (2.24) but can be numerically solved without any difficulty. As in the static case [9], there are expected to be two types of brane solutions, Minkowski and black hole embeddings, depending on whether the brane ends on or above the horizon.
For Minkowski embeddings, the boundary conditions are as before. With these boundary conditions, we find non-trivial profiles f i = f i ( ρ) for i = 0, 4, 8, 10, and f i = 0 otherwise. A general dynamical embedding function is therefore of the form where we neglect terms with i > 11 since they are beyond the validity regime of our approximation of the boost-invariant metric. 3 The numerical plots of the non-trivial profiles f i (i = 0, 4, 8, 10) are shown in figure 3. These profiles have the following qualitative properties. As compared to the solution (2.24) for B = 0, there are additional non-trivial profiles, f 0 and f 4 . The profile f 0 agrees with that in the static case ( τ → ∞) [15][16][17] (dashed curve in figure 3a) We found that the equation for f 0 is independent of ε 0 and η 0 , which is natural since at τ → ∞ the system achieves equilibrium at low temperature and will not depend on the non-equilibrium dynamics (η 0 ) or specific initial conditions ( ε 0 ) anymore. The profiles f 4 and f 8 depend on ε 0 and have a negative sign, which reflects the attraction of the D7 brane to the black hole.  f 4 , f 8 , and f 10 will be suppressed at large τ by negative powers of τ . For self-consistency we will only consider the τ region where all the sub-leading terms are well dominated by the leading terms, i.e. f 0 ≫ f 4 τ −3/4 ≫ · · · . Some of the final embedding profiles L( τ , ρ) are shown in figure 4, where the green lines are plotted by plugging the numerical data of figure 3 into (3.15). We would not expect this expansion approach to work beyond the point where the symmetry breaking embedding ceases to be even a local minimum because higher order terms will grow. One can nevertheless plot solutions using f 0 − f 10 that reach down all the way to the horizon. We will use the full PDE solutions in the next section to test the point where the expansion has broken down.
In fact early in our studies we tried to use the expansion even for black hole embeddings. The ODEs for the f i do not contain the horizon however. We attempted to put the horizon in by hand by imposing boundary conditions relevant to a black hole. For each point in the L − ρ plane if we assume a horizon is present we can deduce the value of τ from (3.11) -one can impose some boundary condition such as orthogonality to the posited horizon and seek solutions for each f i . This correctly gives one massless solution at each τ and they look very similar to the equilibrium black hole embeddings. In fact though on solving the full PDEs we realize that this is far beyond the point where the expansion method has collapsed! We will include some of the resulting curves below though as evidence of the break down of the ODE approximation. Figure 4 shows the embedding profiles at various stages of the evolution of the expanding plasma. The black hole horizon (3.11) is indicated by a black quarter circle in the figure. Its size decreases with time corresponding to a cooling and expanding plasma. Figure 4 reflects the division of the quark-gluon plasma into the three phases described above.
The quark condensate as a function of time can be read form the asymptotic form of JHEP01(2011)050 the solutions We plot these results in figure 2a (the dashed lines) for comparison to the equilibrium inspired results. The main result here is that the late time behaviour is indeed just the equilibrium expectations. Where the curved line deviates from the equilibrium results is in fact a sign that the expansion used in this section has broken down, the configuration is out of equilibrium and adiabatic behaviour is no longer possible. Full solutions of the PDEs will show this in the next section.

Full PDE solutions
In the previous section we interpreted our solution obtained from the ordinary differential system as the evolution of extremum states in the potential. This evolution also reflects the real time dynamics of the plasma whenever there is a unique deep global minimum for the vacuum state, because the embedding time dynamics is expected to be well localized around the minimum. This is the case in the chiral symmetric phase at early times and in the broken phase at late times. In the mixed phase at intermediate times there appears also a local minimum in addition to the global minimum. The small potential barrier separating these minima is potentially easily overcome by a fluctuation. The path of the local maximum seen in the adiabatic approach will not be realized as a dynamical solution since it is unstable. Further below the critical time τ a where the symmetry breaking minimum disappears a heating vacuum will be left in a very non-equilibrium state that again can not be followed using the expansion technique of section 3.3.
To find the real time evolution of the chiral transition out of equilibrium, we need to solve the PDE directly and compare the solution with our previous approximate solutions. In practice it is more convenient to consider a heating process than a cooling one 4 because we can use a well defined starting configuration at large τ as an initial condition of our partial differential equation, i.e.
with f 0 as in (3.15). For simplicity, we impose only Neumann boundary conditions at ρ = 0 associated to Minkowski embeddings and a zero bare quark mass condition at ρ → ∞ to study the spontaneous symmetry breaking,  Numerically we solve the equation using Mathematica's inbuilt PDE solvers. These are somewhat temperamental and one needs to spend considerable time adjusting precision tolerances in order to find smooth solutions in sensible periods of computer time. When we have such solutions we test their stability to changes in precision settings to ensure they are reliable.
With the boundary conditions above we can run our simulations until the embeddings touch the black hole. Beyond that one needs dynamic boundary conditions along the black hole surface. At least in the coordinates we use here this is a hard problem. We have found a simple trick that seems to produce sensible black hole embeddings though. After the D7 has touched the horizon at ρ = 0 we artificially hold the embedding at the top of the horizon. The large ρ evolution of the D7 is local and relatively unaffected by this incorrect embedding at small ρ. Further, as has been observed before, solutions shooting from the black hole horizon experience a numerical attraction onto the unique regular black hole embedding ending at a given point on the horizon. The result is that we get numerical solutions like those shown in figure 5 (e.g. the red, blue, or orange curve) where the D7 follows the horizon before shooting out to large ρ. We believe that these represent very good approximations to the large ρ embeddings solutions. Since we extract the condensate c at large ρ we will live with the improper near horizon behaviour. It would of course be interesting to try to improve on this with dynamic boundary conditions in the future. As a first example of a solution we will study the super-heated symmetry breaking vacuum. At large τ we use the leading terms in the expansion from section 3.2 to find the UV configuration -one needs to use several terms in the expansion to find the embedding with no extra energy. We then solve for the evolution to low τ . In figure 5 we show plots of the embeddingL( ρ) for different τ with the black hole's position for each τ also shown. We expect that the near horizon behaviour is not correct but the far UV embedding should approximate the solution we seek well. In the 3d plot of L( τ , ρ) shown in figure 6 this corresponds to excising the interior of the region indicated by the dashed blue line.
A smooth evolution is apparent. To compare this to our various approximations above we also plot the condensate as a function of time in figure 2a (blue curve). Again the solution follows the equilibrium estimate and the expansion solution at large τ . In the period τ a−c it follows the equilibrium result not the ODE expansion results showing that expansion had broken down before the brane left the mixed phase. The success of the equilibrium approximation suggests we should take its estimate of the transition points τ a (where the local symmetry breaking minimum vanishes) and τ c (where the two minima of the mixed phase are degenerate in energy) as correct. The full PDE solutions allow us to know in addition the behaviour of the condensate when we have heated above the temperature where the super-heated phase has stopped having a local minima (beyond τ a ). This is the main result of our analysis here.
It would be nice to test the equilibrium configurations estimate for the length of time in which the mixed phase exists ( τ b − τ a which is about 10 units of τ ). This we can do by looking at some simple out of equilibrium configurations. For the plots so far shown we computed˙ L( ρ) at some large τ from the f 4 term in the expansion for L( ρ, τ ). We can give the configuration more energy by simply multiplying that˙ L( ρ) by a numerical factor, κ. 5 One might think of these initial states as thermally excited versions of the asymptotic vacuum. In figure 7 we plot the evolution of L( ρ = 0) as a function of τ for a number of such states. At large τ when the theory is cold the unique vacuum is the symmetry 5 In practice, we fix some e τmax and modify the initial conditions breaking one and with extra energy the configuration oscillates about that minimum. The motion is simple harmonic as can be seen from the independence of the period on the amplitude of the oscillations. As the solutions approach τ c the different solutions begin to diverge. The solution with κ = 10 lies close to the equilibrium curve -the oscillations about the minimum are small and the state stays super-heated. When κ = 30, 40 L(0) falls more quickly suggesting that at least some of the brane's length in ρ has escaped the local potential minimum. The upwards wiggles suggest that some of the length is still repelled back into the well though. Finally though by κ = 50 the brane has ridden over the potential barrier and escaped the local minimum. The difference in arrival times at the horizon for these configurations (about 6 units of τ ) is a rough estimate of the period of the mixed phase. It seems to broadly match the equilibrium inspired picture again.
In figure 8 we plot the motion of the IR ρ = 0 end point of the brane with time for a low energy configuration (κ ≃ 1). The brane certainly seems localized in a minimum down to a time of order τ = 27 (compare to the equilibrium estimate that the mixed phase begins at τ b = 34.9 and becomes meta-stable at τ c = 27.5). Further the steepest period of acceleration is below τ = 22 (to be compared with the equilibrium estimate for when the metastable vacuum ceases to exist τ a = 20.6). We conclude that this system is indeed a super-heated state that survives in the local minimum until very close to the equilibrium estimate for τ a . Note the other obvious feature in the plot is the deceleration just below τ = 20. This corresponds to where in our simulation the D7 first impacts on the black hole -at smaller τ we hold the D7 at the horizon as discussed above so this behaviour is an artefact. A full solution would continue to accelerate along the horizon.
These configurations with excessL are also very much linked to bubble formation. A bubble forms in the mixed phase when a volume of space has excess energy due to a thermal fluctuation and escapes the local minima early. Here by treating the whole space as one we are essentially describing the formation of a large homogeneous bubble. It would be interesting in the future to try to study x dependent initial conditions to understand how quickly or slowly bubbles grow. In the previous section we focused on the chiral transition induced by a magnetic field, which simply played the role of an intrinsic symmetry breaking scale. In this section we turn to the physics depending on B more quantitatively.
To make analytic progress we will concentrate on the adiabatic (ODE) approach in which we rescaled all variables by some power of B, see (3.9). In order to study the dependence of the quark condensate on B, we use the original parameters, in terms of which the condensate can be expanded as Note that c 0 is independent of ε 0 , η 0 and η 0 enters only in c 10 . The B dependence of the leading term agrees with that in the static (zero temperature) case [15][16][17], i.e. it scales as B 3/2 . The first subleading term (∼ τ −4/3 ) may be compared to the finite temperature case [15-17, 22, 23]. In the adiabatic approximation, where T ∼ τ −1/3 , this term scales approximately like T 4 B −5/6 . However, due to the B dependence of c 4 , this scaling is not exact. In general, the effect of the subleading terms is to lower the exponent in c(τ, B) ∼ B ν to a value ν < 3/2. The scaling of the total condensate c(τ, B) with B will again be determined numerically.  Figure 10. τ * as a function of B (ε 0 = 1) from the ODE adiabatic approximation results.
in figure 2a). We take τ * to mark the time of the phase transition. As can be seen from figure 2a the early time behaviour is an imperfect approximation to the full PDE solutions but the later time solutions are consistent with the solutions. We find that for fixed τ , or equivalently for fixed temperature, the condensate grows with increasing magnetic field. In the limit τ → ∞ (i.e. at zero temperature), we find the following dependence on B: which is in agreement with the zero temperature result [15][16][17]. For earlier times, the dependence on B is shown in figure 9b. We find numerically that c(τ, B) ∼ B ν where the power ν decreases from 1.5 at large τ to approximately 1.0 at small τ . In other words, the dependence on B tends to become linear at high temperatures. Our results hold for sufficiently strong magnetic fields. At small B the system is in the symmetric phase (c = 0). The tendency for the condensate to increase with B is in qualitative agreement with observations in chiral perturbation theory [38,39] (∝ B 3/2 for strong fields, ∝ B for weak fields), in the Nambu-Jona-Lasinio model [40] (∝ B 2 ), in a confining deformation of the holographic Karch-Katz model [19] (∝ B 2 ), and in SU(2) [41] (∝ B) and SU(3) [42] lattice calculations (∝ B 1.6±0.2 ). The dependence on B typically ranges from linear to quadratic behaviour, i.e. the powers of B are in the range 1 ≤ ν ≤ 2.
We may also study the effect of B on the time of the chiral transition which is marked by τ * . From the D-brane picture we expect that τ * decreases with B. For large B, the repelling force caused by the B-field is much stronger than the attractive force of the black hole. Even at early times, Minkowski embeddings, associated to stable mesons, are therefore favoured over black hole embeddings, which implies that the meson melting process sets in at some earlier time, i.e. at higher temperatures. Figure 10 shows τ * as a function of B. We find that τ * indeed decreases with increasing B as τ * (B) ∼ B −1.55 . The numerically found exponent −1.55 for the scaling of B is close to −1.5 for the case of vanishing shear viscosity, which can be explained as follows.  Figure 11. The condensate c( τ ). The green curve is the leading term, the blue is up to subleading term, the red is up to third term and the black is the whole expression.
At τ = τ * and η = 0, the horizon is located at where the last equality is from the numerical value ofr * at ε 0 = 1. Thus where for the numerical analysis we chose ε 0 = 1. The deviation from −1.5 is due to the effect of the shear viscosity. We also numerically confirmed that the exponent is 1.5 without viscosity. This scaling of τ * may be compared to results for the critical temperature T c of the phase transition in the static case. In the adiabatic approximation (when η 0 = 0), T ∼ τ −1/3 , which implies T * ∼ B 1/2 . This square root behaviour is in agreement with the result for the critical temperature T c in the static approach [22,23] and with [39] for strong magnetic fields. It is also in qualitative agreement with studies in QCD [43,44].
Finally we consider the effect of changing the viscosity. The viscosity effect is very small, since it is doubly suppressed by both large τ and small η 0 . This can also be seen in figure 11, where the condensate is plotted for four cases: the green is the leading term, the blue includes the subleading term, the red includes up to the third term and the black is the full expression. The viscosity effect is then the difference between the red and black curves. There is a small visible difference only at small τ .

Discussion
We have analyzed the first order chiral and meson melting phase transition in a warming or cooling strongly coupled gauge theory with quarks using the AdS/CFT Correspondence. We have developed numerical techniques to study the PDE that describes the motion of a D7 brane in a time dependent geometry. In particular this allows us to explicitly find smooth solutions of the non-equilibrium configurations that are necessarily part of the transition. These results confirm the equilibrium analysis of the transition but also go beyond them. For example we have described the formation of large homogeneous bubbles JHEP01(2011)050 in the mixed phase of the transition. In the future it would be interesting to study spatially inhomogeneous bubbles to understand their growth although this would require the solution of a 2+1 dimensional PDE which is potentially more numerically tricky.
To keep our problem a simpler 1+1 dimensional PDE we also restricted motion of the D7 brane to the holographic coordinates L and ρ. There could also of course be fluctuations in the holographic angular direction φ which we have not described. Such configurations might be useful in the study, for example, of disordered chiral condensates [45]. We may also consider the finite density case by turning on the time component of the U(1) gauge field on the probe brane [22,23,46,47].
We hope that the system we have studied can shed some light on first order transitions in a range of strongly coupled gauge theories including perhaps QCD. The N = 4 background, however, does deconfine in the presence of an infinitesimal temperature, which is not expected in simple QCD. We note though that there are several arguments in favour of a possible splitting of the deconfinement and chiral transitions in QCD in the presence of a strong magnetic field [43]. While the temperature of the chiral transition increases [43,48], the temperature of the deconfining transition either decreases with increasing magnetic field [49] or increases but much slower than the chiral transition [43]. Both transitions become first-order transitions and a new phase with broken chiral symmetry and deconfinement appears for sufficiently strong magnetic fields. So, it is the region of the T -B phase diagram with strong magnetic fields and above the deconfinement temperature in which our model might be qualitatively compared with QCD. (1−3w 4 ) 2 X + +(L ′ ) 2 +1

(A.5)
Open Access. This article is distributed under the terms of the Creative Commons Attribution Noncommercial License which permits any noncommercial use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.