Universal far-from-equilibrium Dynamics of a Holographic Superconductor
MMIT-CTP 4553LA-UR-14-24054
Universal far-from-equilibrium Dynamics of aHolographic Superconductor
Julian Sonner a , † Adolfo del Campo,
2, 3, 4 and Wojciech H. Zurek CTP, Laboratory for Nuclear Science,Massachusetts Institute of Technology, Cambridge,77 Massachusetts Avenue, MA 02139, U.S.A. Department of Physics, University of Massachusetts, Boston, MA 02125, USA Theoretical Division, Los Alamos National Laboratory, Los Alamos, NM 87545, USA Center for Nonlinear Studies, Los AlamosNational Laboratory, Los Alamos, NM 87545, USA
Abstract
Symmetry breaking phase transitions are an example of non-equilibrium processes that requirereal time treatment, a major challenge in strongly coupled systems without long-lived quasiparti-cles. Holographic duality provides such an approach by mapping strongly coupled field theories inD dimensions into weakly coupled quantum gravity in D+1 anti-de Sitter spacetime. Here, we useholographic duality to study formation of topological defects – winding numbers – in the course of asuperconducting transition in a strongly coupled theory in a 1D ring. When the system undergoesthe transition on a given quench time, the condensate builds up with a delay that can be deducedusing the Kibble-Zurek mechanism from the quench time and the universality class of the theory,as determined from the quasinormal mode spectrum of the dual model. Typical winding numbersdeposited in the ring exhibit a universal fractional power law dependence on the quench time, alsopredicted by the Kibble-Zurek Mechanism. a after 1 July 2015: D´epartement de Physique Th´eorique, Universit´e de Gen`eve, 1211 Gen`eve 4, Switzerland † [email protected] a r X i v : . [ h e p - t h ] J un . INTRODUCTION Non-equilibrium quantum phenomena are of wide importance across several disciplinesof physics. Despite their fundamental relevance, few widely applicable principles are knownfor field theories far from equilibrium [1–3]. Gauge-gravity duality is a powerful tool in thisrespect, as it provides a well-defined first-principles framework to study strongly-coupledfield theories in the non-equilibrium setting. These theories are strongly correlated in thesense that there exists no weakly-coupled quasi-particle picture upon which one could basea perturbative treatment. This makes understanding their dynamics, even near equilibrium,an extremely challenging problem. Consequently progress has been confined mostly to casesin which the dynamics is integrable [4–7].Holography [8] is a powerful tool enabling us to explicitly analyze such theories withouthaving to rely on integrability: One maps the quantum-field theory of interest to a dualgravity problem in asymptotically Anti de Sitter (AdS) space-time, which can be solvedin great detail. This duality can reveal significant new insights, for example it providesnew examples of interaction driven localization transitions [9], as well as modelling finite-temperature transport near the superfluid/insulator critical point [10].Previous work on dynamics of holographic superfluids analyzed the condensation pro-cess by perturbing an uncondensed initial state below criticality [11], as well as the non-equilibrium phase diagram resulting from sudden quenches of the superfluid phase [12]. Bothof these studies focus on spatially homogeneous configurations, as do [13, 14], which presentscaling results for the timescale of the breakdown of adiabaticity during finite-rate quenches.Subsequent work, relaxing the constraint of spatial homogeneity, explored vortex turbulenceof holographic superfluids [15] (see also [16] for recent work on anisotropic quenches). Thisarticle extends the framework in yet another direction, namely into the regime of non equi-librium dynamics across phase transitions. We analyze the breaking of a continuous U (1)symmetry, by studying the evolution of black holes in AdS, giving a dual description ofsuperconductors. Here we explore a canonical non-equilibrium paradigm associated withbroken symmetries, leading to universal scaling results in the formation of defects when acritical point is crossed at a finite rate [17, 18]. In this context, inhomogeneous configura-tions of the order parameter are known to be crucial during the formation of topologicaldefects. 2 ote: Independently, Chesler, Liu and Garcia-Garcia [19] have explored the Kibble-Zurek scaling via the dual black hole quasi-normal mode spectrum and numerical analysisin two spatial dimensions.
II. RESULTSA. Winding number generation
Phase transitions were traditionally studied as equilibrium phenomena. Thus, in thebroken symmetry phase, the whole system was assumed to make the same selection of thebroken symmetry vacuum. However, as noted in the cosmological context, where a sequenceof symmetry breaking phase transitions is thought to have resulted in the familiar fundamen-tal forces, rapid cooling of the post Big Bang Universe combined with relativistic causalitymakes it impossible to coordinate symmetry breaking outside of the Hubble horizon. As aconsequence, distinct domains of the Universe will choose broken symmetry vacua on theirown. The resulting mosaic of broken symmetry vacua will – in the course of the subse-quent phase ordering – attempt to smooth out. As Kibble pointed out [17], these disparatechoices can crystalize into topological defects that may have significant consequences for thesubsequent evolution of the Universe.As noted by one of us [18], systems undergoing second order phase transitions cannot everfollow a sequence of instantaneous equilibria. This is because of the critical slowing downin the vicinity of the critical point: the relaxation time τ of the order parameter diverges asa function of the dimensionless distance (cid:15) from e.g., the critical temperature T c : (cid:15) = T − T c T c , (1)where T denotes the instantaneous temperature of the system. This implies that the “re-flexes” of the system are characterized by the universal power-law, τ ( (cid:15) ) ∼ τ | (cid:15) | − zν , (2)where z and ν are the dynamic and correlation length critical exponents, respectively and τ is a microscopic parameter. An arbitrary cooling ramp as a function of time, t , can belinearized around the critical point (cid:15) = − tτ Q . (3)3s a result of critical slowing down, the system loses the ability to adjust to the changeeven when it happens slowly, on any finite quench timescale τ Q . The instant ˆ t when thesystem can no longer keep up with the change of (cid:15) happens when its relaxation time becomescomparable with (cid:15)/ ˙ (cid:15) , the rate of change of (cid:15) . This leads [18] to the equation: τ ( (cid:15) (ˆ t )) = (cid:15)/ ˙ (cid:15) = ˆ t. (4)Using the critical slowing down scaling relation, Eq. (2), one obtains:ˆ t = (cid:0) τ τ νz Q (cid:1) νz . (5)This time scale allows us to split the crossing of the phase transition into a sequence of threestages in which the dynamics is first adiabatic ( t < − ˆ t ), then effectively impulse duringthe interval ( − ˆ t, ˆ t ), and finally adiabatic again deep in the broken-symmetry side of thetransition ( t > ˆ t ). The value of ˆ (cid:15) corresponding to the time scale separating the frozen andadiabatic stages, ˆ (cid:15) = (cid:18) τ τ Q (cid:19) νz , (6)enters into the sonic horizon estimate, the analog of the causal horizon in the early universe.As a result, the characteristic size ˆ ξ of the domains that can coordinate the choice of brokensymmetries exhibits a universal power-law dependence [18] on the quench time:ˆ ξ = ξ ˆ (cid:15) − ν = ξ (cid:18) τ Q τ (cid:19) ν νz , (7)where ξ is a microscopic parameter. The density of topological defects can then be estimatedby recognizing that a ˆ ξ -sized fragment of defect can be expected within a volume ˆ ξ -sizeddomain, the sonic horizon size. Here the velocity of the relevant sound assumes the roleof the speed of light in the relativistic cosmological setting. This reasoning is expected toyield correct scaling of the density of defects, but only an order of magnitude estimate of theprefactor. The scenario just described is often referred to as the ‘Kibble-Zurek mechanism’and we will refer to it throughout this paper as ‘KZM’.Here we consider the setting where the phase transition happens in a ring of circumference C . In view of our above discussion one can expect that the broken symmetry will bechosen independently in sections of size ˆ ξ , so there will be C/ ˆ ξ fragments of the ring thatselect broken symmetry independently [18]. Consequently, phase mismatch resulting from4 random walk of phase with C/ ˆ ξ steps needed to circumnavigate C is expected to leadto: ∆Θ ≈ (cid:113) C ˆ ξ . The net “phase distance” ∆Θ will then have to settle to a multiple of 2 π ,defining the winding number, W , as W = ∆Θ2 π ≈ π (cid:115) C ˆ ξ . (8)It follows that the dispersion of the values of W will scale with the quench rate as [20, 21]: σ ( W ) = 12 (cid:115) C ξ (cid:18) τ Q τ (cid:19) − ν zν ) . (9)This universal power-law is expected to hold whenever C/ ˆ ξ (cid:29)
1, away from the onset ofadiabatic dynamics. The prefactor has to be taken with the proverbial “grain of salt”. Asseen in previous numerical experiments, such KZM calculations yield correct scalings, buttend to overestimate the density of defects [20, 22–24] as well as the typical winding numbers[21, 25].
B. A superconducting ring in holography
We study a one-dimensional superconducting ring of circumference C . The quantity ofinterest is the winding number W = (cid:73) C d Θ( φ )2 π ∈ Z , (10)where φ is the angle along the ring. It follows from gauge invariance that in equilibriumwe have (cid:72) ( d Θ − A ) = 0 so that a non-vanishing winding number is accompanied by a non-vanishing line integral of the vector potential around the loop. It should be noted that W as defined above is gauge invariant under single-valued gauge transformations.We consider the so-called ‘bottom-up’ holographic superconductor model [26, 27] in three-dimensional AdS [28], denoted from here on AdS and work in the probe limit. This amountsto neglecting the effect of the charged components of the system on the neutral ones [15,26]. From the gravity point of view this corresponds to neglecting the bulk gravitationalbackreaction of the charged scalar and the Maxwell field, but keeping the backreaction of thescalar and Maxwell field on each other. This means that we consistently study the dynamicsof a Maxwell field A coupled to a scalar field ψ on a fixed gravitational background. We5 BH bcde <|O(t)|> a rt c FIG. 1.
Schematic of the holographic mapping.
The dual space-time is a cylinder, with timerunning vertically upwards. The field theory lives on the mantle, also called the boundary andshown here in red, an infinite proper distance from the bulk. In the interior resides a black hole(BH) with temperature C/ πz h , which is also the temperature of the dual field theory. Coolingcorresponds to the change of the black hole radius z h . For convenience all equations in this paperuse the coordinate z ∈ (0 , z h ], which is related to the more conventional AdS radial coordinatevia (cid:37) = L /z ∈ [ (cid:37) h , ∞ ). The radial coordinate r is a compactified version of (cid:37) , which puts theboundary at a finite distance. Panel a schematically shows the order parameter as function oftime, spatially averaged over C . The four time slices, b - e , indicate the condensation process asit happens throughout the bulk. The density profile shown in each slice schematically illustratesthe behavior of the bulk field ψ dual to the order parameter in the boundary theory with blueindicating a vanishing density of | ψ | and red a high density of | ψ | . b : At the beginning thereis no condensate, but the temperature is starting to be lowered through T c , which is reachedat time t c , as indicated in a . c : Sometime after crossing T c a small, spatially inhomogeneouscondensate starts appearing. d : The condensate grows, amplifying the original inhomogeneitiesinto macroscopic phase domains. e : Eventually the condensate settles down to its equilibriumconfiguration with a given winding number frozen in. choose the scalar field to have vanishing mass, which means that it is dual to a classicallymarginal operator in the boundary field theory. We do not expect our results to differ6ualitatively for other choices of the mass, as long as it remains in the range correspondingto a marginal or relevant operator (we refer here to the renormalization group, ‘RG’, scalingdimension in the ultra-violet, ‘UV’). The details of our bulk action and equations of motioncan be found in the Methods section. Maxwell theory without the charged scalar in AdS andits dual field theory have been studied previously in [29–31], taking advantage of the fact thatone may conveniently dualize the bulk vector field to a bulk scalar. A different holographicapproach to superconducting rings, using probe branes, was developed in [32, 33].Studying the field theory at finite temperature corresponds to studying the bulk system inthe background of a black hole, in our case the three-dimensional black hole (the BTZ blackhole, after Ba˜nados, Teitelboim and Zanelli, [34]), which is characterized by the parameter z h , its horizon radius. For more details see Fig. 1. With this data we have the Hawkingtemperature T H = 12 πz h . (11)Via the AdS/CFT correspondence (CFT stands for conformal field theory) this is directlytranslated into the temperature of the dual field theory [35]. This theory is therefore studiedin a thermal ensemble at temperature T H , which is an external parameter in our dynamics.In order to fully specify the ensemble (in the equilibrium case), as well as the dynamicswe must augment the equations of motion with suitable boundary conditions at the UVboundary z = 0. We impose in addition to the finite temperature T a fixed charge density ρ corresponding to Neumann boundary conditions on the bulk gauge field. We give Dirichletboundary conditions for the scalar field. From this it follows that the field theory is studiedin the absence of a source for the order parameter density. This condition ensures thatany symmetry breaking occurring will be spontaneous. In order to compare our full non-equilibrium results to the prediction of KZM, we first need to determine the universalityclass of the field theory just defined, as well as the microscopic parameters τ and ξ . As isusual in holography, this boils down to computing the so-called quasinormal modes of thebulk black hole [36]. C. Near equilibrium universality
In order to determine the critical exponents z, ν , as well as the microscopic parameters τ and ξ , it is sufficient to study the holographic superconductor near equilibrium. We will7 b ct /τ Q = .
72 t /τ Q = .
82 t /τ Q = . Θ Θ Θ r (cid:45) r (cid:45) r (cid:45) r (cid:45) r (cid:45) r (cid:45) d e fg h i | ψ | | ψ | | ψ | B B B
FIG. 2.
Winding up a superconductor following a temperature quench.
Example con-densation process with τ Q = 12 leading to a W = 3 configuration, shown in three stages. Wechoose a gauge that makes the final phase of the order parameter linear Θ( φ ) = W φ . The top row(panels a , b , c ) shows the field-theory condensation process; the local magnitude of the condensateis represented as the radius of the torus, while the local phase is encoded as the color, rangingfrom [0 : red] to [2 π : blue]. The middle row ( d - f ) shows the magnitude | ψ | of the bulk fielddual to the order parameter O while the last row ( g - i ) shows the magnetic field B in the bulk.The bulk magnetic fields disappear at late times and the winding around the boundary circle isaccompanied by winding number around the horizon of equal magnitude. study the non-compact case, since the finite size of the ring in our simulations does notaffect the results for local correlations, so long as the healing length ξ is much smaller than8he circumference. The non-compact spatial boundary direction is denoted x here. Since ξ formally diverges near T c , this assumption will break down for extremely slow quenches, butwe have not found that our simulations entered this regime for the range of τ Q under study.Adapting the analysis of [37] to the present situation, the correlation function of the orderparameter field O near, but slightly above, T c takes the momentum-space form (cid:104)O ( ω, k ) O † ( − ω, − k ) (cid:105) := χ ( ω, k ) = Z ( ω, k ) icω + k + ξ , (12)where the last expression holds for small ω, k and we introduced the parameters c and 1 /ξ .We have defined the dynamical susceptibility of the order parameter χ ( ω, k ) in the firstequality. One obtains the correlation function from the susceptibility χ ( ω, k ) by a Fouriertransform G ( t − t (cid:48) , x − x (cid:48) ) = (cid:90) dωdk (2 π ) e iω ( t − t (cid:48) ) − ik ( x − x (cid:48) ) χ ( ω, k ) . (13)We can find the relaxation time by looking at the Fourier transform of the zero-momentumresponse, χ ( ω, k = 0), which at late times takes the form G ( t − t (cid:48) ) ∝ e − t − t (cid:48) cξ , t (cid:29) t (cid:48) , ⇒ τ = cξ . (14)The equal time correlation function, following from the Fourier transform of the static cor-relation function χ ( ω = 0 , k ), falls off like G ( t = t (cid:48) , x − x (cid:48) ) ∝ e − x − x (cid:48) ξ , x (cid:29) x (cid:48) , ⇒ z = 2 . (15)From the fact that the relaxation time is proportional to the square of ξ it follows that thedynamical critical exponent z = 2. This leaves us to determine the correlation length criticalexponent ν . For this we must study how the correlation length ξ diverges near the criticalpoint, ξ ∼ ξ | (cid:15) | − ν . Using the above relations (12) & (14), we deduce that τ ∼ τ | (cid:15) | − ν with τ = cξ . In holography the poles of two-point functions of boundary operators coincide withthe quasinormal modes of the bulk fields dual to the operators appearing in the correlationfunction [38]. As we saw above, in order to extract ξ and τ it is sufficient to study thestatic susceptibility χ ( ω = 0 , k ) and the dynamical susceptibility χ ( ω, k = 0) separately. Asshown in Fig. 3, we find that the spectrum of quasinormal modes of χ ( ω, k = 0) containspoles only in the lower half complex plane (as demanded by stability), while χ ( ω = 0 , k )has two series of poles along the positive and negative imaginary axis. This structure is also9pparent from the correlation function (12). In each case the relevant poles are the onesclosest to the real axis, governing the exponential decay at long time scales (defining therelaxation time τ ), and the falloff of correlations at large distances (defining the correlationlength ξ ), respectively. We find from the quasinormal mode analysis that 1 /ξ ∼ | (cid:15) | , whichimplies ν = 1 /
2, consistent with the results of [37]. By studying the motion of the leadingpoles as T → T c as discussed in Methods, we can also deduce that τ = 2 . ± .
01 and ξ = 0 . ± .
02 if the critical point is approached from above and τ = 0 . ± .
01 and ξ = 0 . ± .
01 if it is approached from below. For more details, we refer the reader toFig. 3. 10 ⊙⊙⊙⊙⊙⊙⊙⊙⊙⊙ ×× - ��� - ��� ��� ��� ��� - � - ���� �� ( � / ρ ) � � ( � / ρ ) b ���� ���� ���� ���� ���� ���� ���� ���������� � / � � ρ ξ c üü üü ¥¥ - - - - - - - - - H w ê r L I m H w ê r L d ���� ���� ���� ���� ���� ��������������� � / � � ρ τ �� FIG. 3.
Susceptibilities in the normal phase.
Panels a & c : Poles in correlation functionsat T /T c = 1 .
1, as deduced from studying the quasinormal modes of the bulk black hole. Thesemodes are linear fluctuations with dissipative boundary conditions at the horizon. The leadingpoles nearest the real axis, marked by red crosses, determine the behavior of correlations at largedistances ( a ) and long times ( c ). Momentum k and frequency ω are given in units of charge density ρ , while temperature T is given in units of the critical temperature T c . Panels b & d : divergenceof correlation length ( b ) and critical slowing down ( d ), as determined from the correspondingpoles of χ ( ω, k ) shown in a and c . The dashed curve in b shows a fit to ξ = ξ | (cid:15) | − ν resulting in ξ = 0 . ± .
02 and ν = 0 . ± . d shows a fit to τ = τ | (cid:15) | − zν , resultingin τ = 2 . ± .
01. Furthermore, we can analytically determine z = 2. Within numerical accuracy,we find identical values for z and ν below the critical point, with slightly different ξ = 0 . ± . τ = 0 .
89. We have the approximate relation ξ < = √ ξ > between ξ below and above thetransition. This is the opposite of the usual Landau-Ginzburg mean-field relation ξ < = 1 / √ ξ > .Thus, while the critical exponents have the mean-field values, equilibrium correlations lengthsclearly show that the theory we are dealing with is outside of the the Landau-Ginzburg paradigm. . Far from equilibrium dynamics We simulate cooling of the superconducting ring by numerically evolving the bulk equa-tions (Eqs. (20) in Methods) in the black-hole background with changing temperature (11).We implement a piece-wise linear protocol starting at an initial temperature of T i with cor-responding (cid:15) i , and cooling the system at finite rate according to (cid:15) ( t ) = − t/τ Q through thecritical point, until the temperature T f with corresponding (cid:15) f is reached. We implement thetemperature ramp by changing the dimensionless ratio T C = C/ πz h , while holding the den-sity ρC fixed. For a precise definition of these quantities the reader may consult the Methodssection. In order to allow the system to break the symmetry dynamically, one needs to addfluctuations to the classical Einstein equations. To achieve this, we introduce noise into theevolution in a manner consistent with the fluctuation-dissipation theorem for the horizontemperature T H . Thus we sample the boundary values, ψ ( z, t, φ ) (cid:12)(cid:12) z = z UV , of the scalar field(its real and imaginary parts) from a Wiener process, which ensures that their average valuesvanish but its dynamics gives rise to a non-vanishing correlation. We treat the amplitudeof the noise as a phenomenological parameter, but it would be enlightening to determineits precise form in the future, for example by deriving the relevant fluctuation-dissipationrelation from bulk quantum effects (for work in this direction, see [39–41]). Before the lin-early decreasing temperature ramp we allow the system to thermalize for a fixed time inorder to allow the noise introduced at the boundary to get distributed over all scales in thebulk. In order to solve the nonlinear partial differential equations determining the evolutionof the bulk fields, we use a Chebyshev grid in the holographic (‘radial’) direction z and aFourier decomposition in the periodic direction φ . We use a fixed step to evolve in time.Each quench is started at some (cid:15) > ˆ (cid:15) for a given τ Q and is stopped at a time, t stop whenthe condensate, averaged over the ring, has reached a fixed fraction of the equilibrium valueat (cid:15) ( t stop ). The time elapsed between the time of crossing the equilibrium phase transitionpoint, and the stopping time t stop defines the lag time t L . At this point the winding number W is recorded together with the lag time t L . For the simulations of this work, a fraction of0 . W is extracted from the discretized φ direction (denote the grid points12y { φ i } N x i =1 ) via the sum W = 12 π N x (cid:88) i =1 Arg (cid:16) O i +1 O † i (cid:17) . (16)As in the case of the continuum definition above, the sum as a whole is gauge invariantunder single-valued gauge transformations. We find that it stabilizes soon after the orderedphase is formed (see Fig. 4). Winding number becomes a good observable by the time westop the evolution to extract both t L and W , when it is frozen in at integer values. Fromthis point on it is no longer susceptible to noise or late-time equilibration dynamics of theordered phase, as can be seen in Fig. 4. - ��� - ��� ��� ��� ��������������������� � / τ � < � ( � ) > / < � ( ∞ ) > - ��� - ��� ��� ��� ��� ��������������������� � / τ � < � ( � ) > / < � ( ∞ ) > - ��� - ��� ��� ��� ��� ��������������������� � / τ � < � ( � ) > / < � ( ∞ ) > - ��� - ��� ��� ��� ��� - � - � - ����� � / τ � � � ∫ � / ( � π ) - ��� - ��� ��� ��� ��� ��� - � - � - ����� � / τ � � � ∫ � / ( � π ) - ��� - ��� ��� ��� ��� ��� - � - � - ����� � / τ � � � ∫ � / ( � π ) a b cd e f FIG. 4.
Time-evolution of the order parameter and the winding number.
Representativesamples of individual runs for three different quench rates. First column (panels a & d ): τ Q = 10;second column (panels b & e ), τ Q = 35; third column (panels c & f ): τ Q = 120. In the first row ( a - c ) we show the time development of the order parameter (red curve) compared to the instantaneousadiabatic value (green-black dashed). In the second row ( d - f ) we show the time evolution of thewinding number (black) compared to the line integral of the gauge field π (cid:72) A . The two convergein equilibrium, as dictated by gauge invariance: (cid:72) d Θ = (cid:72) A . The winding number is always integerquantized, the non-vertical parts of the curve are an artifact of joining up the quantized values atfinite sampling intervals. The winding number is physically well defined only once a condensatehas developed. An example condensation process for τ Q = 12 is shown in Fig. 2, where boundary andbulk physics leading to a W = 3 configuration is illustrated.13e find very good agreement between the full simulations and KZM predictions based onthe equilibrium critical exponents deduced from our independent quasinormal mode analysis.From ν = 1 / z = 2 it follows that (cid:104)| W |(cid:105) ∝ τ − / and t L ∝ τ / , while the simulationresults in (cid:104)| W |(cid:105) ∝ τ − . ± . , t L ∝ τ . ± . . (17)The quoted uncertainties give a single standard deviation from the fitted value. The dis-crepancy in the value of t L is likely a result of the fact that the lag time does not vanish atvery small quench times, but rather saturates to a finite value, so that the simple scalingform is no longer a very good fit for the rapid quenches at the fast end of our window of τ Q .In summary, we find good agreement with universal KZM values for the scaling of t L as wellas the dispersion of winding number σ ( W ). Matching the prefactors in (7) with equilibriumpredictions usually results in more significant quantitative departures. In the present casewe obtain ξ sim0 = 1 . ± . , τ sim0 = 3 . ± . , (18)deviating by a factor of ∼ ξ = 0 . τ = 2 .
02) extractedfrom correlation functions. From past experience it is to be expected that the KZM valuesover-estimate the number of defects. KZM predictions of σ ( W ) for τ Q = 10 , , . , . , .
37 (rounded to two significant digits), compared to 1 . , . , .
00 (again,rounded) from the full simulation, so KZM overestimates the density of defects. Evidentlythe numbers in our simulation are in rather good agreement, compared to past simulationswhere mismatches by factors of O (10) were not uncommon (see e.g. [21]). Provided thatthe numerics is good enough, the degree of agreement depends on the microscopic dynamics.Furthermore our value for τ is in line with the mismatch encountered in previous work [21].14 τ Q σ ( W ) b (cid:1) (cid:2)(cid:3) (cid:4)(cid:3) (cid:1)(cid:3) (cid:2)(cid:3)(cid:3) (cid:1)(cid:3)(cid:3) (cid:2)(cid:3)(cid:3)(cid:3)(cid:1)(cid:2)(cid:3)(cid:1)(cid:3)(cid:2)(cid:3)(cid:3) (cid:3)(cid:5)(cid:3) (cid:3)(cid:5)(cid:6) (cid:3)(cid:5)(cid:7) (cid:3)(cid:5)(cid:8) (cid:3)(cid:5)(cid:9) (cid:2)(cid:5)(cid:3) (cid:2)(cid:5)(cid:6)(cid:3)(cid:5)(cid:3)(cid:3)(cid:5)(cid:2)(cid:3)(cid:5)(cid:6)(cid:3)(cid:5)(cid:4)(cid:3)(cid:5)(cid:7) τ Q t L c ��� ��� ��� ��� ��� ��� ��������������� t/ ˆ t (cid:104) O ( t ) (cid:105) FIG. 5.
Universal scaling laws for the lag time and the dispersion of winding number.
Panel a : Best fit to winding number scaling gives excellent agreement with values predicted byKZM. Shown is a fit to σ ( W ), which gives σ ( W ) = aτ bQ , with a = 2 . ± .
25 and b = − . ± . γ τ − / , where equilibrium estimates for γ give thevalue 3 .
99. The error bars in panel a were computed as the standard deviation of the mean foreach τ Q in our dataset for winding numbers. The error bars in panel b give a single standarddeviation from the mean for each τ Q in the dataset. We see that the dispersion in winding numbersaturates for rapid quenches. These deviations are accompanied by a deviation from the naivescaling prediction of ˆ t at the onset of saturation for σ ( W ) (see also panel c and its descriptionbelow). Panel b : Best fit to scaling of the lag time t L . Inset shows the condensate as a functionof t/ ˆ t for the three cases of Fig. 4. Panel c : Order parameter O ( t ) averaged over the ring as afunction of time, given in units of ˆ t . This last quantity is computed from the near-equilibriumcritical behavior of our system. We show four different values of the quench time, τ Q = 10 in black, τ Q = 35 in red, τ Q = 120 in blue, and finally τ Q = 300 in green. The scaling collapse for lag timebecomes more and more accurate for slower quenches, in accordance with the results in panels a & b , showing more significant deviations for fasts quenches. II. DISCUSSION
The work reported here constitutes the first demonstration of universal KZ scaling ofdefects spontaneously created in the far-from equilibrium dynamics of a strongly-correlatedfield theory without quasiparticles. Holography allowed us to map the time-dependenceof this system to the evolution of a set of nonlinear partial differential equations togetherwith stochastically sampled boundary conditions. Our setup is comparable to a holographicversion of the Gross-Pitaevskii equation, applicable for systems without quasiparticles.We introduced stochasticity by sampling boundary conditions from a Wiener processcharacterized by a phenomenological parameter α . This is natural, because it constitutesonly a minimal extension of the usual holographic dictionary, which tells us how to mapbulk to boundary quantities, but has the limitation that there is a free parameter. While wefound that the universal results concerning the scaling law as a function of the cooling rateare very robust to changes in α , and the parameters τ sim0 and ξ sim0 of our study show onlyweak dependence, it would nevertheless be desirable to derive it from first principles. Thismeans deriving the fluctuation-dissipation relation associated with the thermal Hawkingradiation of the bulk black hole, which will force us to take into account quantum effects inthe bulk gravity model along the lines of [39–41].Our scaling results are consistent with exponents obtained in the mean-field approxima-tion of the boundary theory, even though the theory is not described by a simple Landau-Ginzburg (LG) model, as evidenced e.g. by a violation of the standard LG relation betweencorrelation length above and below the transition (see Fig. 3 and its caption). The mean-field like scaling results are a consequence of working in the classical gravity limit, which,however, does manage to capture some features beyond mean-field, notably the ratio of cor-relation lengths above and beyond the phase transition. In holography, this quantity takesthe form of a ratio of bulk integrals [42], and thus receives contributions from all scales,including the horizon, the source of dissipation in our holographic representation. The factthat dissipative dynamics are naturally incorporated into the theory from first principlesmakes holography a powerful tool to investigate strongly-coupled superfluid and supercon-ductor dynamics, as previously stressed in [15], where the holographic origin of dissipationwas crucial in establishing a turbulent direct cascade. In more conventional approaches,such as the stochastic Gross-Pitaevskii or Landau-Ginzburg approaches, dissipation has to16e added by hand. A tempting next goal would be to investigate quench dynamics in the-ories that do not result in mean-field scalings. A much more ambitious (but perhaps notcompletely unrealistic [43]) goal would be to study strongly coupled theories directly relevantto condensed matter physics.An open question left unanswered by the present investigation is the origin of the satu-ration we observed. There are several plausible explanations, and we hope to return to adetailed investigation of the reasons for the saturation in the model we have presented herein forthcoming work.Another interesting future project would be to study the effect of the noise in more detail.In addition to the first principles derivation of the properties of the Wiener process we havealready noted, two distinct phenomenological tacks can be considered. We have alreadybegun to investigate the effect of the amplitude of the Wiener process (our parameter α )applied at the boundary. The preliminary conclusion is that – at least within the range wehave explored – the quantities we have followed are insensitive to suitably small α , with thedependence likely to be sublogarithmic. More detailed characterization of this dependencewould be desirable. Furthermore, one can consider noise applied throughout the AdS interior(rather than only on the boundary). This raises the possibility of seeding topological defectsin the interior, and may pose the question of the relation between them and the behavior ofthe field on the boundary. Acknowledgements:
We would like to thank Allan Adams, Tarek Anous, Chris Herzog, Nabil Iqbal, KristanJensen, Arttu Rajantie, Homer Reid and Toby Wiseman for very helpful discussions. Thenumerical computations in this paper were performed on the MIT LNS openstack cluster,and we thank Jan Balewski and Paul Acosta for their kind assistance. We also thank AllanAdams for providing us with further computational resources. This research is supportedby the U.S Department of Energy through the LANL/LDRD Program and a LANL J.Robert Oppenheimer fellowship (AdC). This work was also supported in part by the U.S.Department of Energy (DOE) under cooperative research agreement Contract Number DE-FG02-05ER41360 (JS). 17
V. METHODSA. Bulk action and Details of the Holographic Mapping
The full action for the AdS holographic superconductor reads S = 116 πG N (cid:90) d x √− g (cid:18) R + 2 (cid:96) (cid:19) − q (cid:90) √− gd x (cid:18) F + | Dψ | + m | ψ | (cid:19) , (19)where D µ ψ = ∇ µ ψ − iA µ ψ and we choose the case m = 0. We denote bulk indices runningover z, t, x by x µ and boundary indices, running over t, x , with x i . Taking q → ∞ has theresult of decoupling the gravity part from the matter part [26], tantamount to the ‘probelimit’ discussed in this article. One fixes a solution of the gravity equations, to be discussedshortly, and treats the dynamics of the matter fields separately. Hence the effective bulkaction from which the equations of motion follow is simply Maxwell-scalar theory in a non-trivial background, that is, the equations of motion are D ψ = 0 , ∇ µ F µν + J ν = 0 , (20)for the current J µ = i ( ψ ( D µ ψ ) ∗ − ψ ∗ D µ ψ ). Covariant derivatives are taken with respect toa fixed background metric, which we take to be the metric of the BTZ black hole, writtenin ingoing Eddington-Finkelstein coordinates ds = (cid:96) z (cid:0) − f ( z ) dt − dtdz + dx (cid:1) . (21)The coordinate x is identified x ∼ x + C for the dynamical simulations, where for numericalconvenience a rescaled coordinate φ = 2 πx/C is used. The Eddington-Finkelstein timecoordinate t reduces to boundary time at z = 0, warranting the use of the same symbol forboth. The coordinate x is kept unidentified for the calculations of the quasinormal modesin section D. The length (cid:96) sets the AdS curvature scale; its role is to fix the number ofdegrees of freedom ‘ N ’ of the dual theory, which must be large for classical gravity to apply.At leading order in large N this quantity scales out of the equations we solve. The metricfunction takes the fixed form f = 1 − z /z . For the scalar field one finds the asymptoticexpansion ψ ( t, z, x ) = ψ (0) ( t, x ) + ψ (2) ( t, x ) z + · · · (22)18he interpretation of ψ (0) ( t, x ) is that it sources the symmetry-breaking operator O ( t, x ),while ψ (2) ( t, x ) gives its expectation value ψ (2) = 12 (cid:96) (cid:104)O(cid:105) . (23)Therefore the requirement that the symmetry be broken spontaneously means that we mustset the source ψ (0) ( t, x ) to zero for all time. This translates into a homogenous Dirichletboundary condition on the field ψ ( t, z, x ) in the UV, ψ ( t, z, x ) (cid:12)(cid:12)(cid:12) z = z UV = 0 . (24)The vector field in AdS is more subtle [29]. Its asymptotic behavior is given by A µ ( t, z, x ) = j µ ( t, x ) log( z/ Λ) + a µ ( t, x ) + · · · , (25)where the vector j µ is an external current in the boundary theory, not to be confused with thebulk current J µ ( t, z, x ). We introduced the scale Λ to make the argument of the logarithmdimensionless. We shall see that our chosen boundary condition on A µ is independent ofthis scale. Note that in the normal phase the solution is j µ log( z/z h ) and Λ = z h is enforcedby regularity at the horizon. It is convenient to work in axial gauge so that A z = 0 andthus the current j µ ( t, x ) has components only in the field theory directions, j i ( t, x ). Notethat this is not the same as choosing axial gauge in the Schwarzschild like coordinate system(39), in which the bulk metric is diagonal. In the latter choice of coordinates, the equationsof motion imply the equation ∂ i j iS − i(cid:96) (cid:0) ψ (0) (cid:104)O(cid:105) ∗ − ψ ∗ (0) (cid:104)O(cid:105) (cid:1) = 0 . (26)Thus the current is conserved in the absence of a source for the operator O ( t, x ) (We havedenoted the current evaluated in the Schwarzschild like coordinates as j iS ). From the point ofview of the boundary field theory this is simply the Ward identity for the one-point functionof the current. Operationally, the conservation condition follows from the z − component ofthe Maxwell equations, which gives rise to a constraint. Turning to the gauge field, theboundary condition is Π µA (cid:12)(cid:12)(cid:12) z = z UV + j µ = 0 ⇒ z(cid:96) ∂ z A µ (cid:12)(cid:12)(cid:12) z = z UV = j µ , (27)19here Π µA is the momentum conjugate to A µ with respect to z slicingΠ µA = lim z → z UV √− gF µz . (28)This is a Neumann boundary condition. Thus we are free to fix j i ( t, x ) subject to theconservation condition, and leave a i ( t, x ) free to fluctuate. As stated above, the scale Λ dropsout from the boundary condition (27). This choice corresponds to a dual vector operator ofdimension ∆ = 1, the right dimension for a dynamical gauge field in the boundary theory.Indeed the residual gauge transformation preserving axial gauge ( λ res ( t, z, x ) = λ ( t, x )) actson this as a standard field-theory gauge transformation a i ( t, x ) → a i ( t, x ) + ∂ i λ ( t, x ) , (29)so that a i ( t, x ) is a bona-fide fluctuating gauge field. We fix a constant background chargedensity, so that the current has the only non-vanishing component, j t = ρ . B. Details on Numerics
The simulations in this article were performed on a pseudo-spectral spatial grid comprisedof 21 Chebyshev points in the radial direction and 111 plane waves in the angular directionof the boundary.For each run we start the system at the initial time slice in the normal phase, defined bysetting to zero the field ψ and giving the gauge field a non-trivial time component j t = ρ .We then evolve forward in time, updating the noise according to the rule (16), implementedas a discrete Wiener process. We average over O (10 ) noise realizations for each value of τ Q to compute σ ( W ). Computer codes used in the simulations of this work are available uponrequest. C. Evolution Scheme
For the simulations reported on in this article we employed a characteristic evolutionscheme for the Maxwell and scalar fields. After gauge fixing A z = 0, denoting A t ( t, z, x ) = T ( t, z, x ), A x ( t, z, x ) = X ( t, a, x ) and writing ψ ( t, z, x ) = a ( t, z, x ) + ib ( t, z, x ) we havefour evolution equations and one constraint equation. In the numerical evolution it proves20onvenient to work with rescaled fields, rather than the ‘bare’ ones appearing in the action(19) . We furthermore subtract the leading log terms. Thus defining ˆ T = z ( T + ρ log z )and ˆ X = zX we find the equations (cid:20) ∂ z − z (cid:21) Φ it = S i [ a, b, T, X ] , (cid:20) ∂ z − z (cid:21) ˆ T t = S T [ a, b, T, X ] , (30)where the fields Φ i are all fields other than ˆ T . The sources S i and S T depend non-linearlyon the fields as well as their spatial derivatives. In addition, the radial ( z ) component of theMaxwell equation gives the constraint equation (cid:0) − z∂ z + z ∂ z (cid:1) ˆ T = 2 (cid:96) z ( b∂ z a − a∂ z b ) + z ∂ zx ˆ X − z∂ x ˆ X . (31)After the rescaling, the boundary conditions on ˆ T and ˆ X are nowˆ T ( t, z, x ) (cid:12)(cid:12)(cid:12) z =0 = ˆ X ( t, z, x ) (cid:12)(cid:12)(cid:12) z =0 = 0 . (32)Note that Eqs. (30) are linear equations for the time derivatives suggesting the followingevolution scheme. Assume that at time t = t n the values of the fields a, b, ˆ X are known.We can now use the constraint equation to solve for ˆ T ( t n , z, x ) and then solve the linearequations for x it to obtain the time derivatives of a, b, ˆ X . We can then use our favorite timeevolution scheme to obtain the values of a, b, ˆ X on the next time slice t n +1 . We have usedexplicit RK4 integration as well as simple forward Euler with good results. For 111 Fouriermodes in the spatial boundary direction, the step size was taken to be 0 . /N z , where N z is the number of grid points in the radial direction. The scaling analysis of Figs. 5 wasobtained on a grid of 21 points in the radial direction and 111 Fourier modes in the annulardirection choosing a characteristic size of the ring C = 50 (cid:96) . It would be desirable to repeatthe analysis for a larger ring with higher spatial resolution, which is likely to require moresignificant computing power, especially if higher statistics on noise realizations is desired.The linear equations for ∂ t Φ in on a given time slice can then be solved efficiently andin parallel for each spatial grid point φ j . For the above-mentioned spatial grid we choosea fixed time step of ∆ t = 0 . /N z . In order to achieve stable long-time evolution wefilter out high-frequency modes in the φ direction using the Orszag 2 / D. Holographic Renormalization
In axial gauge, and setting ψ (0) = 0, the equations of motion (20) have solutions withasymptotic behavior A i ( z ; t, x ) = j (0) i ( t, x ) log ( z ) + a (0) i ( t, x ) + · · · ψ ( z ; , t, x ) = z ψ (2) ( t, x ) + · · · (33)We regularize the on-shell action by introducing a finite cutoff z = (cid:15) , so that we now have S = − (cid:90) z = (cid:15) d x √− γ ˆ n µ (cid:2) g µν (cid:0) ( D ν ψ ) ψ † + ψD ν ψ † (cid:1) + F µν A ν (cid:3) + E . O . M . + · · · , = − (cid:96) (cid:90) z = (cid:15) dtdx (cid:104) j (0) i ( t, x ) log (cid:15) + j (0) i ( t, x ) a (0) i ( t, x ) (cid:105) + · · · (34)where ˆ n is the outward pointing unit normal to the boundary and γ = det( γ ij ) the de-terminant of the induced metric. The omitted terms depend only on the boundary values j i , a i , ψ (2) and vanish in the limit (cid:15) →
0. The divergent terms can be cancelled by adding acounterterm S CT = 12 (cid:96) (cid:90) d x √− γF zi F zi log (cid:15) . (35)As shown in [30] this counterterm becomes a manifestly local boundary quantity whenwritten in terms of a dual scalar field S , obtained from the two-form gauge field F via F = (cid:63)dS . E. Noise
Under classical evolution, in the case at hand corresponding to the classical large N limit of the dual field theory, fluctuations are suppressed and the symmetry of the order22arameter cannot be dynamically broken. Said differently, even though below T c the phasewith (cid:104)O(cid:105) = 0 is unstable, there is nothing in the classical evolution equations to push theorder parameter off its precarious perch on top of the potential. Taking our lead fromthe literature on BEC dynamics using the stochastic Gross-Pitaevksii equation [20] we addfluctuations by sampling our boundary conditions from a thermal noise distribution, that isfrom a Wiener process. Thus the boundary conditions for a ( t, z, x ) and b ( t, z, x ) should bemodified. Instead of imposing the strict Dirichlet boundary condition (24) , we only imposethis condition on average, i.e. we impose (cid:104)(cid:104) a (0) ( t, x ) (cid:105)(cid:105) = (cid:104)(cid:104) b (0) ( t, x ) (cid:105)(cid:105) = 0, with (cid:104)(cid:104) a (0) ( t, x ) a (0) ( t (cid:48) , x (cid:48) ) (cid:105)(cid:105) = α TT c δ ( t − t (cid:48) ) δ ( x − x (cid:48) ) , (cid:104)(cid:104) b (0) ( t, x ) b (0) ( t (cid:48) , x (cid:48) ) (cid:105)(cid:105) = α TT c δ ( t − t (cid:48) ) δ ( x − x (cid:48) ) , (36)where (cid:104)(cid:104)·(cid:105)(cid:105) denotes noise average and ψ (0) = a (0) + ib (0) . Usually one determines α from afluctuation-dissipation relation as α = 2 ηT c , where η is a damping parameter. In this workwe treat α as a phenomenological parameter. We found that varying α , even by several ordersof magnitude has little influence on the scaling results, but we do find a weak dependenceof the absolute magnitude of t L on α . Clearly it is desirable in future to determine α froma first-principles holographic calculation.In practice the noisy boundary condition is realized by sampling each spatial boundarypoint from an independent normal distribution of zero mean and unit variance N (0 , a (0) ( t n , x i ) = α TT c √ ∆ tN (0 , ,b (0) ( t n , x i ) = α TT c √ ∆ tN (0 , , (37)for each boundary point x i . Note that one cannot choose both the boundary value of j i and ψ (0) independently, since they are constrained by the current Ward identity (26) . Thisequation is consistent only for ψ (0) = 0, so in what sense can one choose the noisy boundarycondition above? This can be seen by considering the noise as a small perturbation j i + δj i and ψ (0) + δψ (0) satisfying ∂ i δj i − i(cid:96) (cid:0) δψ (0) (cid:104)O(cid:105) ∗ − δψ ∗ (0) (cid:104)O(cid:105) (cid:1) = 0 . (38)By solving the z component of the Maxwell equation near the boundary we automaticallyimpose this constraint at each time step, that is the noise fluctuations (37) also imply afluctuation in the current δj i , such that (38) is satisfied.23ote in particular that (37) implies that the phase is a random variable. In our numericswe find that it is sufficient to update the noise boundary condition at larger intervals, sayevery 100 time steps. We have also run the simulations updating the noise at every timestep, as well as every 10 time steps, again with no apparent impact on the results - providedthe noise amplitude is adjusted in accordance with (37). Clearly there is a lower limit on thesampling frequency for which this statement is true - consider, e.g. , the extreme case of onlyupdating the noise once or twice during the entire simulation. However our results show nodetectable dependence on sampling frequency, which implies that we stayed far away fromthis lower limit throughout. An illustration of the dependence of the evolution on the noiseparameter is given in Supplementary Fig. 1. F. Quasinormal Mode Analysis
We find it convenient to perform this analysis using Schwarzschild coordinates, in whichthe background metric takes the form ds = (cid:96) z (cid:20) − f ( z ) dτ + dz f ( z ) + dx (cid:21) , (39)with the same function f = 1 − z /z h used throughout and Schwarzschild time τ . Wethen expand the full equations of motion to linear order around a background solution. Animportant simplification is that in order to determine poles in the correlation function ofthe order parameter we only need to consider scalar operators, since these cannot mix withvector modes to first order in perturbation theory.We start with an analysis of the normal phase. Infinitesimal gauge transformations of thebackground configuration induce first-order changes in the metric scalar and gauge fields.In order to extract physical poles of correlation functions it is most convenient to determinethe minimal set of gauge-invariant fluctuations containing the mode of interest (here theorder-parameter fluctuation). In the broken-symmetry phase this requires more work, butin the symmetric phase the analysis is very simple. A suitable set of gauge-invariant modesare given by the real and the imaginary part of the scalar field fluctuation around zerobackground values a ( z, τ, x ) = e − iωτ + ikx α ( z ) , b ( z, τ, x ) = e − iωτ + ikx β ( z ) . (40)24he modes α and β decouple from all other fluctuations. This decoupling even happenswhen one goes beyond the probe limit and allows fluctuations of the metric, see, e.g. [12].A more detailed discussion of the gauge-invariant modes for the probe system is presentedbelow in the context of the broken phase. The modes we just identified satisfy the equations α (cid:48)(cid:48) + (cid:18) f (cid:48) f − z (cid:19) α (cid:48) + 1 f (cid:0) ω − f k + T (cid:1) α + 2 iωT f β = 0 ,β (cid:48)(cid:48) + (cid:18) f (cid:48) f − z (cid:19) β (cid:48) + 1 f (cid:0) ω − f k + T (cid:1) β − iωT f α = 0 . (41)Since the equations are linear in ω , but quadratic in k , these admit an expansion [37] insmall frequency and momentum α = α (0 , ( z ) + ωα (1 , ( z ) + k α (0 , ( z ) + · · · ,β = β (0 , ( z ) + ωβ (1 , ( z ) + k β (0 , ( z ) + · · · , (42)where α ( i,j ) ∼ α s ( i,j ) + z α v ( i,j ) + · · · β ( i,j ) ∼ β s ( i,j ) + z β v ( i,j ) , (43)for z → z UV . In this expression a superscript s denotes source behavior, while v denotesexpectation value. Thus Green functions of the dual operator have the small ω, k behavior G ( ω, k ) = α v (0 , + ωα v (1 , + k α v (0 , α s (0 , + ωα s (1 , + k α s (0 , := Z ( ω, k ) icω + k + 1 /ξ , (44)as claimed in (12) above. More precisely, the Green functions of the operators dual to α and β mix, but the eigenvalues of the matrix of correlation functions will be of the functionalform (44). This establishes analytically that the dynamical critical exponent is z = 2. Inorder to determine ξ , and thus ξ and ν , we solve the system of equations (41) numerically,using a finite difference discretization. Results of this analysis are described in more detailin Sec. II. above, and summarized in Fig. 3.We now describe the more involved calculation in the broken phase, where for ω (cid:54) = 0and k (cid:54) = 0 more mode mixing occurs. If we were to go beyond the probe approximation,the sound channel [44], holographically encoded as a scalar fluctuation of the metric, willalso contribute. This means that in general the condensate-condensate two points functionbelow T c no longer takes the simple form (12).Below the critical temperature, i.e., when the background contains a nontrivial scalarfield, we have to take into account the modes e − iωτ + ikx { a τ ( z ) , a z ( z ) , a x ( z ) , α ( z ) , β ( z ) } , where25 µ ( z ) are perturbations of the gauge field and α ( z ) and β ( z ) are perturbations of thereal and imaginary parts of the complex scalar field respectively. Not all of these modesare physical, since we can generate an infinitesimal perturbation by acting on the back-ground with a gauge transformation e − iωτ + ikx λ ( z ). This generates the perturbation modes { δa τ , δa x , δa z } = e − iωτ + ikx {− iωλ ( z ) , − ikλ ( z ) , λ (cid:48) ( z ) } in the gauge sector, as well as a pertur-bation of the imaginary part of the scalar δβ = e − iωτ + ikx a ( z ) λ . We cannot generate a realpart of the scalar perturbation in this way, because the background value a ( z ) is purelyreal. Thus we have three gauge-invariant perturbationsΦ = α , Φ = iωβ + a a τ , Φ = − ikβ + a a x . (45)The last mode, Φ transforms as a vector, but the combination k Φ is a scalar, so that it toocan contribute to the linear equations for the condensate fluctuation. Since this combinationvanishes at zero momentum, it does not appear in the dynamic susceptibility calculation. Dynamic Susceptibility: k = 0 In this case the vector like perturbation Φ decouples, and we have the equationsΦ (cid:48)(cid:48) + (cid:18) f (cid:48) f − z (cid:19) Φ (cid:48) + W [ T ]Φ + W [ T , a , ω ]Φ = 0 , Φ (cid:48)(cid:48) + V [ T , a , ω ]Φ (cid:48) + W [ T ]Φ + W [ T , a , ω ]Φ = 0 . (46)with the coefficients W = ω + T f , W = 2 iωT f , W = − iT ωz f (cid:0) ω z − a f (cid:1) ,V = ω z ( zf (cid:48) − f ) − a f ( a − za (cid:48) ) zf ( ω z − a f ) ,W = 1 z f ( ω z − a f ) (cid:104) a f + ω z (cid:0) ω + T (cid:1) − z a f (cid:0) ω + T (cid:1) − z f a (cid:48) − za a (cid:48) f ( zf (cid:48) − f ) (cid:105) . (47)For completeness we also present the decoupled equation for the remaining gauge-invariantmode Φ (cid:48)(cid:48) + (cid:18) z + f (cid:48) f (cid:19) Φ (cid:48) + ω z − a fz f Φ = 0 . (48)The behavior of these poles in the complex plane is summarized in Supplementary Fig. 2.26 tatic Susceptibility: ω = 0 For the computation of the static susceptibility, Φ decouples from the other two modes,which satisfy the coupled equationsΦ (cid:48)(cid:48) + (cid:18) f (cid:48) f − z (cid:19) Φ (cid:48) − k f − T f Φ + 2 a T f Φ = 0 , Φ (cid:48)(cid:48) + 1 z Φ (cid:48) − k z + 2 a z f Φ − a T z f Φ = 0 . (49)The remaining equation for the mode Φ readsΦ (cid:48)(cid:48) + (cid:18) f (cid:48) f − z + 4 a ( a − za (cid:48) ) z ( k z + 2 a ) (cid:19) Φ (cid:48) + (cid:18) a (cid:48) ( za (cid:48) − a ) z ( k z + 2 a ) + T f − z k + 2 a z f (cid:19) Φ = 0 . (50)The behavior of these poles in the complex plane is summarized in Supplementary Fig. 3. Author Contribution:
The authors jointly defined the project. JS developed and carried out the numerical sim-ulations. All authors contributed to the analysis and interpretation of the numerical dataand the preparation of the manuscript.
Additional Information
Competing Financial Interests:
The authors declare no competing financial interests.27 ppendix A: Supplementary Material a ��� ��� ��� ����� - � �� - � �� - � ������������� (cid:104) O ( t ) (cid:105) / (cid:104) O ( ∞ ) (cid:105) t/τ Q b ���� ���� ���� ���� ���� ���� ���� ������ × �� - � �� × �� - � �� × �� - � �� × �� - � ������� (cid:104) O ( t ) (cid:105) / (cid:104) O ( ∞ ) (cid:105) t/τ Q Supplementary FIG. 1.
Noise dependence of the order parameter dynamics.
Panel a:Average over the ring of condensate density as a function of time for two different values of thenoise amplitude differing by a factor of ten, shown in red and black. Panel b shows a blowup ofthe ‘knee’ region where the condensation process first happens. We see that this process is largelyinsensitive to the noise amplitude. üü ü ¥¥ - - - - - - - - H w ê r L I m H w ê r L b ê T c t eq Supplementary FIG. 2.
Dynamic susceptibility in the broken phase.
Panel a : Poles indynamic correlation function at T /T c = 0 .
9. There is a Goldstone pole at the origin. The leadingpole nearest the real axis, corresponding to the amplitude or ‘Higgs’ mode, with imaginary part ω (cid:63) gives the equilibration time. Panel b : inverse of the imaginary part of leading pole ω (cid:63) as thecritical point is approached. This directly gives the equilibration time τ . The best fit result gives τ ∼ . (cid:15) − , with a τ differing from the unbroken phase. Again the explicit calculation agrees withthe analytical derivation that z = 2. The pole structure found here resembles closely the AdS results of [12, 45]. üüüüüüüüü ¥¥ - - - - - H k ê r L I m H k ê r L b ê T c x Supplementary FIG. 3.
Static susceptibility in broken phase.
Panel a : poles in staticcorrelation function at T /T c = 0 .
9. The leading poles nearest the real axis, with imaginary part k (cid:63) give the correlation length. Panel b : imaginary part of leading pole k (cid:63) as the critical point isapproached from below. The correlation length ξ = k − (cid:63) diverges as ξ = ξ (cid:15) − / . We determined ξ = 0 . ± .
01 from numerically solving Eq. (41).[1] Dziarmaga J. Dynamics of a quantum phase transition and relaxation to a steady state,
Adv.Phys. , 1063–1189 (2010).[2] Polkovnikov, A., Sengupta, K., Silva, A. & Vengalattore, M. Colloquium: Nonequilibriumdynamics of closed interacting quantum systems, Rev. Mod. Phys. , 863–883 (2011).[3] del Campo, A. & Zurek W. H. Universality of phase transition dynamics: Topological defectsfrom symmetry breaking. Int. J. Mod. Phys. A , 1430018 (2014).[4] Calabrese, P. & Cardy, J. Phys. Rev. Lett. , 136801 (2006).[5] Calabrese, P., Essler, F. H. L. & Fagotti, M. Quantum quench in the transverse-field Isingchain. Phys. Rev. Lett. , 227203 (2011).[6] Fioretto, D. & Mussardo, G. Quantum quenches in integrable field theories.
New J. Phys. ,055015 (2010).[7] Cuax, J.-S. & Konik, R. M. Constructing the generalized Gibbs ensemble after a quantumquench. Phys. Rev. Lett. , 175301 (2012).
8] Maldacena, J. M. The Large N limit of superconformal field theories and supergravity.
Adv.Theor. Math. Phys. , 231–252 (1998).[9] Donos, A. & Hartnoll, S. A. Interaction-driven localization in holography. Nature Phys. ,649–655 (2013).[10] Witczak-Krempa, W., Sørensen, E. S. & Sachdev, S. The dynamics of quantum criticality viaquantum Monte Carlo and holography. Nature Phys. , 361–366 (2014).[11] Murata, K., Kinoshita, S. & Tanahashi, N. Holographic superfluids and the dynamics ofsymmetry breaking. JHEP , 050 (2010).[12] Bhaseen, M. J., Gauntlett, J. P., Simons, B. D., Sonner, J. & Wiseman, T. Non-equilibriumcondensation process in a holographic superconductor.
Phys. Rev. Lett. , 015301 (2013).[13] Basu, P., Das, D., Das, S. R. & Nishioka, T. Quantum quench across a zero temperatureholographic superfluid transition.
JHEP , 146 (2013).[14] Basu, P. & Das, S. R. Quantum quench across a holographic critical point.
JHEP , 103(2012).[15] Adams, A., Chesler, P. M. & Hong, L. Holographic Vortex Liquids and Superfluid Turbulence.
Science , 368-372 (2013).[16] Bai, X., Lee, B.-H., Park, M. & Sunly, K. Dynamical condensation in a holographic supercon-ductor model with anisotropy.
JHEP , 054 (2014).[17] Kibble, T. W. B. Topology of cosmic domains and strings.
J. Phys. A: Math. Gen. , 1387–1398 (1976).[18] Zurek, Z. H. Cosmological experiments in superfluid helium? Nature , 505–508 (1985).[19] P. M. Chesler, A. M. Garcia-Garcia and H. Liu, Defect Formation beyond Kibble-Zurek Mech-anism and Holography. Phys. Rev. X (2015) 2, 021015[20] Sabbatini, J., Zurek, W. H. & Davis, M. J. Phase separation and pattern formation in a binaryBose-Einstein condensate. Phys. Rev. Lett. , 230402 (2011).[21] Das, A., Sabbatini, J. & Zurek, W. H. Winding up superfluid in a torus via Bose Einsteincondensation.
Sci. Rep. , 352 (2011).[22] Laguna, P. & Zurek, W. H. Density of kinks after a quench: When symmetry breaks, how bigare the pieces? Phys. Rev. Lett. , 2519-2522 (1997).[23] Antunes, N. D., Bettencourt, L. M. A. & Zurek, W. H. Vortex string formation in a 3-D U(1)temperature quench. Phys. Rev. Lett. , 2824-2827 (1999).
24] De Chiara, G., del Campo, A., Morigi, G., Plenio, M. B. & Retzker, A. Spontaneous nucleationof structural defects in inhomogeneous ion chains.
New J. Phys. , 115003 (2010).[25] Nigmatullin, R., del Campo, A., De Chiara, G., Morigi, G., Plenio M. B. & Retzker, A.Formation of helical ion chains. Preprint at arXiv:1112.1305 [quant-ph] (2011).[26] Hartnoll, S. A., Herzog, C. P. & Horowitz, G. T. Building a holographic superconductor. Phys.Rev. Lett. , 031601 (2008).[27] Hartnoll, S. A., Herzog, C. P. & Horowitz, G. T. Holographic superconductors.
JHEP ,015 (2008).[28] Ren, J. One-dimensional holographic superconductor from AdS /CFT correspondence. JHEP , 055 (2010).[29] Marolf, D. & Ross, S. F. Boundary conditions and new dualities: vector fields in AdS/CFT.
JHEP , 085 (2006).[30] Jensen, K. Chiral anomalies and AdS/CMT in two dimensions.
JHEP , 109 (2011).[31] Faulkner, T. & Iqbal, N. Friedel oscillations and horizon charge in 1D holographic liquids
JHEP , 060 (2013).[32] Khlebnikov, S. Winding branes and persistent currents.
JHEP , 105 (2013).[33] Khlebnikov, S. Critical current of a superconducting wire via gauge/gravity duality.
Phys.Lett. B , 271 (2012).[34] Ba˜nados, M., Teitelboim, C. & Zanelli, J. The Black hole in three-dimensional space-time.
Phys. Rev. Lett. , 1849-1851 (1992).[35] Witten, E. Anti-de Sitter space, thermal phase transition, and confinement in gauge theories. Adv. Theor. Math. Phys. , 505–532 (1998).[36] Horowitz, G. T. & Hubeny, V. E., Quasinormal modes of AdS black holes and the approachto thermal equilibrium. Phys. Rev. D , 024027 (2000).[37] Maeda, K., Natsuume, M. & Okamura, T., Universality class of holographic superconductors. Phys. Rev. D , 126004 (2009).[38] Kovtun, P. K. & Starinets, A. O. Quasinormal modes and holography. Phys. Rev. D ,086009 (2005).[39] de Boer, J., Hubeny, V. E., Rangamani, M. & Shigemori, M. Brownian motion in AdS/CFT. JHEP , 094 (2009).[40] Son, D. T. & Teaney, D. Thermal noise and stochastic strings in AdS/CFT.
JHEP , 021 Phys. Rev. Lett. , 091601 (2012).[42] Maeda, K. & Okamura, T., Characteristic length of an AdS/CFT superconductor.
Phys. Rev.D , 106006 (2008).[43] McGreevy, J., Duality with a view toward many-body physics. Advances in High EnergyPhysics , 723105 (2010).[44] Policastro, G., Son, D. T. & Starinets, A. O. From AdS / CFT correspondence to hydrody-namics. II. Sound waves. Duality with a view toward many-body physics.
JHEP , 054(2002).[45] Amado, I., Kaminski, M. & Landsteiner, K. Hydrodynamics of holographic superconductors.
JHEP , 021 (2009)., 021 (2009).