Topological Defects Formation with Momentum Dissipation

We employ holographic techniques to explore the effects of momentum dissipation on the formation of topological defects during the critical dynamics of a strongly coupled superconductor after a linear quench of temperature. The gravity dual is the dRGT massive gravity in which the conservation of momentum in the boundary field theory is broken by the presence of a bulk graviton mass. From the scaling relations of defects number and"freeze-out"time to the quench rate for various graviton masses, we demonstrate that the momentum dissipation induced by graviton mass has little effect on the scaling laws compared to the Kibble-Zurek mechanism. Inspired from Pippard's formula in condensed matter, we propose an analytic relation between the coherence length and the graviton mass, which agrees well with the numerical results from the quasi-normal modes analysis. As a result, the coherence length decreases with respect to the graviton mass, which indicates that the momentum dissipation will augment the number of topological defects.


Introduction
Non-equilibrium dynamics of strongly coupled system is one of the most challenging tasks in theoretical physics. It is not only interesting in cosmology and particle physics, but also interesting within condensed matter physics [1][2][3]. In recent years, the central question in understanding the strongly coupled theory was how the critical dynamics across a phase transition react to a time dependent coupling, i.e., a thermal or quantum quench. We investigate the critical dynamics of a strongly coupled superconducting phase transition within the AdS/CFT correspondence from string theory [4], focusing on the formation and evolution of topological defects after a thermal quench. Spontaneous topological defect formation in such transition is a fundamental component of the celebrated Kibble-Zurek mechanism (KZM) [5][6][7].
The basic idea of KZM is that: As a continuous symmetry is spontaneously broken during a second order phase transition, the evolution of the system stops being adiabatic and enters the "freeze-out" region as a result of the critical slowing down nearby the critical point. Topological defects will form at the interfaces where different symmetrybreaking domains meet and satisfy the "geodesic rule" [8]. KZM predicts a universal power law between the topological defects number and the quench rate. Such mechanism was first introduced by Kibble [5,6] into Cosmology that as a result of relativistic causality, topological defects such as cosmic strings, monopoles, and etc., can form in the subsequent evolution of the Universe. It was later introduced by Zurek into condensed matter physics that the vortex lines may occur in a superfluid when the system is quenched, which is accessible in the laboratory [7]. This scenario was supported by various numerical studies [9][10][11][12][13][14][15], and experiments carried out in a variety of systems [8,[16][17][18][19][20][21][22][23][24][25].
The theory of KZM was proposed in a spatially homogenous background in which the momentum is conserved. However, real materials, which are composed of electrons, atom lattices, impurities and so on, do not possess such kind of momentum conservations (or have conservation only modulo the reciprocal lattice vectors) due to the spatial inhomogeneities [26]. It will be quite important to investigate the KZM in the strongly coupled systems with momentum dissipation in the framework of holography. From holography, there have been a number of models built to study momentum dissipation by implementing the effects of translational symmetry breaking. Several approaches among these are to include the lattice effects [27], or study the charge transport in the background of a spatially modulated bulk [28][29][30][31], or in the presence of impurities directly [32,33]. However, due to the involved technics in such models, another conceptually distinct approach to holographic momentum dissipation was suggested in [34]. This approach does not make use of any conventional mechanism of translational symmetry breaking, but instead provides an effective bulk description of a theory without momentum conservation from massive gravity. Applications of gauge/gravity duality on non-equilibrium dynamics can be found in [35][36][37][38][39][40][41][42]. Previous holographic work on KZM were carried out in [43][44][45][46][47][48][49], without any momentum dissipation in the background.
Holographic massive gravity exhibits momentum dissipation [50] by giving a mass to the graviton. This mass breaks the diffeomorphism invariance of the gravitation theory, via the holographic dictionary, it in turn violates the conservation of energy-momentum in the dual field theory. Many clues have shown that the effect of massive graviton in the bulk gravity can be considered as the effect from the lattice in the dual field theory on the boundary [51,52]. Recently, a nonlinear massive gravity theory has been proposed by de Rham et.al. (dRGT theory) [53,54], and later it is found to be ghost-free [55,56]. Note that, the dRGT massive gravity is thought to be the only healthy theory in Poincaré invariant setups, and there have been numerous studies in dRGT massive gravity which we will not list all of them [57][58][59][60][61].
Our goal in this paper is to explore the effect of the momentum dissipation (equivalently the effect of graviton mass) on the critical dynamics across a second-order phase transition after a linear quench of temperature from the holographic (3 + 1)-dimensional dRGT massive gravity. Quantized superconductor vortices will turn out due to KZM after quench. We will see that the momentum dissipation has little effect on the KZM scaling laws of the vortex number and "freeze-out" time to the quench rate. However, it is interesting to see that the coherence lengthξ of the order parameter at the "freeze-out" time decreases with respect to the graviton mass m from the quasi-normal modes (QNMs) analysis. Physically, it implies that the momentum dissipation will reduce the coherence length. Inspired from Pippard's formula for the coherence length and mean free path [62], we propose an analytic expression between the coherence lengthξ and the graviton mass m. This analytic relation is in good agreement with the results from QNMs analysis. Moreover, we also investigate the relations between the number of vortices n and the graviton mass m from both numerical simulations and analytic approach, and find that they fit very well as m is not too large. The resulting vortex number will increase as graviton mass grows, which in turn indicates that momentum dissipation will enlarge the number of ovrtices.
The outline of the paper is as follows. In Sec.2, we briefly review the KZM, and build up the holographic model in the (3+1)-dimensional dRGT massive gravity in in Sec.3; Sec.4 contains the main analytical and numerical results of our paper to study the generations of vortices and the relations between graviton mass; The conclusions are drawn in Sec.5; In the Appendix A we briefly introduce the dRGT massive gravity.

A brief review of Kibble-Zurek mechanism
KZM is a paradigmatic theory to describe spontaneous generation of topological defects in non-equilibrium critical dynamics [5][6][7]. In the vicinity of critical point of a second-order phase transition, critical opalescence and critical slowing-down implies that the coherence length ξ and the relaxation time τ diverge as in which ξ 0 and τ 0 are constant coefficients depending on the microscopic physics, while ν and z are the static and dynamic critical exponents in equilibrium states, and is the dimensionless distance to the critical temperature, ≡ 1 − T /T c . KZM assumes a linear quench across the critical point as (t) = t/τ Q , in which τ Q is the quench rate (or quench time). The "freeze-out" time occurs at the instantt when the relaxation time of the system τ equals the time scale of the quench to the critical point, then the dynamics approximately freezes (from nearly adiabatic to approximately impulse behavior within the time interval [−t,t]). Therefore, the "freeze-out" time can be readily obtained aŝ The key insight of the KZM is that the average sizeξ of the domains in the broken symmetry phase is set by the the equilibrium coherence length ξ evaluated at the freeze-out timet, i.e.,ξ = ξ(t). Thus, the size of the symmetry-breaking domain is given by, Topological defects will occur at the interfaces between various symmetry-breaking domains if they satisfy the "geodesic rules" [8]. Therefore, the above estimate of theξ can be recast as an estimate for the resulting number of topological defects, n ≈ Ap where A is the volume size of the system, p is the average probability for the domains to form defects and d is the spatial dimension of the system. This relation is verified to be universal in various dimensions and for different kinds of topological defects [8][9][10][11][12][13][14][15][16][17][18][19][20][21][22][23]. In particular, the generated topological defects are vortices for superfluids or superconductors in d = 2 system. We will investigate the KZM for vortices in a strongly coupled superconductor system with background momentum dissipation by virtue of the AdS/CFT technique in the following.

Holographic setup
The gravity background we adopt is the dRGT massive gravity in (3 + 1)-dimensional bulk spacetime [53,54] (Please consult Appendix A for a brief review of dRGT massive gravity).
In order to study the dynamics in this background, we work under the Vaidya-AdS solution with flat space k = 0 and zero charge q = 0 (We have scaled the AdS radius L ≡ 1.).
in which M is the mass parameter of the black hole, c i (i = 1, 2) are constants and m stands for the graviton mass. Location of the AdS boundary is z = 0 while the horizon is z = z + . Hawking temperature of the black brane thus is We work in the probe limit and adopt the conventional Maxwell-complex scalar action for holographic superconductors [63] is the covariant derivative and m ψ is the mass of the complex scalar field Ψ. Thus, the equations of motion can be obtained It is convenient to work in axial gauge by setting A z = 0, and take the ansatz for other fields as Ψ = Ψ(t, z, x, y), A t,x,y = A t,x,y (t, z, x, y).

Holographic renormalization and boundary conditions
In order to solve the above equations, we need to impose suitable boundary conditions at the horizon (z = z + ) and at the AdS boundary (z = 0). Near the boundary z → 0, fields can be expanded as, with ∆ ± = 3 ± 9 + 4m 2 ψ /2 are the conformal dimensions of dual scalar operators on boundary. Following the AdS/CFT dictionary, the interpretation of Ψ 0 is that it sources the symmetry breaking operatorÔ, while a t and a i are the chemical potential and superfluid velocity in the boundary, respectively. The physical meaning of Ψ 1 and b µ will be clear after the holographic renormalization.
From holographic renormalization [64], varying the renormalized on-shell action with respect to the source terms one can achieve the corresponding conjugate variables. In order to obtain finite on-shell action, we should add the counter terms of the scalar fields S ct = d 3 x √ −γΨ * Ψ into the divergent on-shell action, where γ is determinant of the reduced metric on the boundary. We set Neumann boundary conditions for gauge fields in order to get the dynamical gauge fields in the boundary. In order to have a well-defined variation, the surface term S sf = d 3 x √ −γn µ F µν A ν should also be added, where n µ is the normal vector perpendicular to the boundary. Hence, we can get the finite renormalized on-shell action S ren . The expectation value of the order parameter can be obtained from the holographic renormalization by varying S ren with respect to Ψ 0 , resulting in Ô = Ψ 1 .
We impose Ψ 0 = 0 all time in order to have a spontaneous symmetry breaking of the order parameter. This translates into a homogenous Dirichlet boundary condition for the field Ψ in the AdS boundary, i.e., Ψ| z=0 = 0. Expanding the z-component of the Maxwell equations near boundary, we reach ∂ t b t +∂ i J i = 0, which is exactly a conservation equation of the charge density ρ and current J i on the boundary, since from the variation of S ren one can get b t = −ρ and At the horizon, we demand the regularity of the fields. Since the metric component g tt is zero at the horizon, the component A t should be vanishing, while other fields are finite at the horizon.

Numerical schemes
We choose m 2 ψ = −2 without loss of generality and set c 1 = 1, c 2 = −1/2 in order to make the background thermodynamically stable [52,57]. For convenience of the numerical calculation, we set the black hole mass parameter as M = 1 and the location of horizon at z + = 1. From holographic superconductor [63], increasing the chemical potential or the charge density equals decreasing the temperature of the system. We choose to increase the charge density ρ in our paper. In addition, the temperature of the black hole T has the mass dimension one while the charge density ρ has mass dimension two from the dimensional analysis. Therefore, in order to linearize the temperature near the critical point according to KZM, we quench the charge density ρ as where ρ c is the critical charge density for the static and homogeneous holographic superconducting system. We should note that different graviton masses m correspond to different ρ c , which is shown in We quench the system from normal state with initial temperature T i = 1.4T c to superconducting state with the final equilibrium temperature T f = 0.8T c . The size of the boundary system is (x, y) = (50, 50). We take advantage of the Chebyshev pseudospectral method in the radial direction z and use the Fourier decomposition in the (x, y)directions since the periodic boundary condition is imposed along (x, y)-directions. The grids in z-direction is 21 while the grids along (x, y)-directions are 201 × 201. We thermalize the system by adding random seeds for the fields at the initial time. The random seeds are sampled in the bulk by satisfying the white noise distributions s(t, x) = 0 and s(t, x)s(t , x ) = ζδ(t − t )δ( x − x ), with the amplitude ζ ≈ 10 −3 . The system evolves by using the 4th order Runge-Kutta method with time step ∆t = 0.01. Filtering of the high momentum modes are implemented following the "2/3's rule" that the uppermost one third Fourier modes are removed [65].

KZM scalings of vortex number and "freeze-out" time
In this subsection, we will examine the KZM scaling laws in Eqs.(2.2) and (2.4), depending on the graviton mass. We count the vortex number n as the order parameter Ô saturates to its equilibrium value. Scalings between the number n and quench time τ Q under different graviton mass m are exhibited in the left panel of Fig.2. The shaded regions correspond to the standard deviations from the mean values. In the fast quench regime (small τ Q ) the vortex number approximately saturates due to the finite size effect, which is consistent with previous results in condensed matter or holography [43][44][45][46]. However, in the slow quench regime (large τ Q ), the vortex number will decrease with respect to quench time satisfying the scaling relations as n m=0 ≈ 305 × τ −0.504 Q , n m=1.2 ≈ 390 × τ −0.5059 Q and n m=2 ≈ 805 × τ −0.518 Q respectively. From this figure, we can conclude that as m grows, or equivalently momentum dissipation gets stronger, the number of vortices will increase as well. In other words, from the relation n ∝ 1/ξ 2 we can deduce that the momentum dissipation will reduce the coherence length. We will study this relation between coherence length and the graviton mass in the next subsection in details.
The "freeze-out" timet can be reflected by studying the lag time t L that defined as order parameter begins to grow rapidly [14,43,45]. In numerics we operationally set t L as Ô ∼ 0.1 following [14,43,45]. On the right panel of Fig.2 we exhibit the relation between t L and τ Q . The error bars are not shown since they are very tiny. We see that for fast quench the lag time is almost constant for different m's. However, for slow quench, one can read that t L ≈ 6.45 × τ 0.5102 Q for m = 0, t L ≈ 4.347 × τ 0.5014 Q for m = 1.2 and t L ≈ 3.98 × τ 0.4931 Q for m = 2. Therefore, from the two KZM scaling relations in Eqs.(2.2) and (2.4), one can readily evaluate the dynamic critical exponent z and the static critical exponent ν on the boundary as (ν ≈ 0.514496, z ≈ 2.0246), (ν ≈ 0.50732, z ≈ 1.98221) and (ν ≈ 0.510949, z ≈ 1.90386) for m = (0, 1.2, 2) respectively, in which the exponent matches the mean-field theory values with z = 2 and ν = 1/2. This means graviton mass in dRGT massive gravity does not alter the KZM scaling laws, which in turn indicates that the existence of graviton mass does not affect critical exponents, z and ν, in the boundary field theory. The boundary field theory remains as a mean-field theory, which is consistent with the assumption of AdS/CFT correspondence that the boundary is a mean-field theory in large N limit [66].

Quasi-normal modes analysis
From previous section 2, we see that the average size of the symmetry-breaking domainŝ ξ is approximately the coherence length of the order parameter at the "freeze-out" timê t [7]. Fortunately, from holography one can compute the coherence length from QNMs of the scalar field in the bulk by choosing suitable boundary conditions [67,68]. The QNMs correspond to the poles of the retarded Green's function of the dual field theory. One can read off the coherence length ξ from the correlation function (4.1) where k and ω are the momentum and frequency of perturbation modes respectively, whilec is a parameter. QNMs are computed from small fluctuations of the background fields, such that Ψ(t, z, x) = ψ 0 (z) + ψ(z)e −iωt+ikx , where ψ 0 (z) is the static background of scalar field while ψ(z)e −iωt+ikx is small fluctuations with respect to the background. For simplicity, we have assumed the fluctuations as plane waves in the only x-direction. We calculate the QNMs in the symmetric phase, i.e., in the normal phase that the background ψ 0 (z) is vanishing. This is because from KZM, dynamics ceases to be adiabatic and enters an impulse stage within the time interval [−t,t] as we learned from previous section 2. The starting of the symmetry breaking will occur at the instant t ≈ −t, which is in the normal phase before the critical point. This greatly simplifies the complexity of the calculations.
It is convenient to work out the QNMs in the Schwarzschild-AdS metric of the spacetime [43,45]. Substituting the scalar field Ψ(t, z, x) into the EoMs (3.5), (3.6), and extracting the first order of the perturbations ψ(z), we arrive at the decoupled linear equation as in which A t is the background gauge field. We impose the ingoing boundary conditions at the horizon z + = 1, where T is the Hawking temperature of black hole. Vanishing boundary conditions are imposed at the boundary z = 0. We adopt the determinant method [45] to calculate the QNMs and fix the quench rate as τ Q = 20 (other quench rates will have similar results). More specifically, we employ the Chebyshev differential matrix to convert the differentiation with respect to spacetime coordinate into a matrix. Zeros of the determinant of the resulting matrix are the QNMs after imposing the above boundary conditions. The lowest mode in the QNMs makes the most contribution to the perturbations. In practice, the lowest mode is the mode that is closest to the real axis. We can evaluate the coherence length from Eq.(4.1) by setting ω = 0, and arrive atξ = 1/|Im(k * )| where k * is the lowest mode in those QNMs. After solving the Eq.(4.2) one can get a series of QNMs in the complex k-plane, as is shown in the left panel of Fig.3 for m = 2. The red dots, symmetric and closest to the real axis, stand for the lowest poles k * . Other sub-leading poles (in blue squares), all locate in the imaginary axis, are also symmetric to the real axis. Different graviton masses correspond to different series of QNMs, leading to a relation between coherence lengths and m, i.e.,ξ =ξ(m), which are shown as red dots in the right panel of Fig.3.

Analytic relation betweenξ and m
The analytic relation between the coherence lengthξ and the graviton mass m was seldom discussed. In [52], the authors found that the coherence length would decrease with respect to the graviton mass in the system of a holographic Josephson junction. However, they did not propose an analytic relation betweenξ and m. In this subsection, we will seek this analytic relation betweenξ and m and show that this relation matches well with the numerical results from QNMs analysis.
From the discussions in [34,50,51], the graviton mass m effectively plays a role of momentum dissipation in the boundary field theory. This momentum dissipation might be from scatterings of electrons to lattices or impurities, which lead to the assumption in Drude model that the relation between the mean free path l and scattering time τ s is τ s = l/V F , in which V F is the Fermi velocity [26]. In condensed matter physics, the coherence length ξ can be related to the mean free path l from Pippard's formula [26,62], where ξ p is the coherence length of a pure material without any scatterings, and α is a numerical constant. Then from the relation τ s = l/V F in Drude model, we get Fortunately, from holography [61] τ s can be related to the graviton mass m as τ −1 s = γm 2 where γ is a coefficient independent of m. 1 Therefor, we finally reach where a = ξ p and b = ξ p γ/(αV F ). The fitting of this relation (4.6) is exhibited as blue line in the right panel of Fig.3, with fitting parameters as a ≈ 0.7466 and b ≈ 0.1496. From Fig.3 we see that the numerical results (red dots) from QNMs analysis are fitted very well with the analytic relation Eq.(4.6). Therefore, the coherence lengthξ at the "freeze-out" time will decrease as m increases. This also means stronger momentum dissipation in the boundary field will reduce the average size of the symmetry-breaking domains, which in turn implies that the vortex number n would grow with respect to graviton mass according to n ∝ 1/ξ 2 in Eq.(2.4).

Vortex number and graviton mass
We count the vortex number at the final equilibrium state of the system, where the order parameter is well developed. From the above analysis, we expect that the vortex number n would increase with the graviton mass, and satisfying the following relation via Eq.  can see that the values of the coefficient before m 2 are close to each other within 5.6151% errors. Estimates of parameter b from these two methods, i.e., one is from the fitting of relation n ∼ m, the other is from the fitting ofξ ∼ m, match each other very well! This also indicates that our conjecture of the relation between coherence length and graviton mass, viz. Eq.(4.6) is very trustable.
The overall factor in the front of Eq.(4.7) can even provide us the estimate of the probability p, which is the probability for the symmetry-breaking domains to form vortices, in Eq.(2.4). We already knew a ≈ 0.7466 in Eq.(4.6) from the right panel of Fig.3, and the size of the system is A = 50×50 = 2500, thus we can readily get p ≈ 17.82×0.7466 2 /2500 ≈ 0.003973 and n ≈ Ap/ξ 2 ≈ 9.9331ξ −2 . The later relation is plotted as red dashed line in the Fig.4, whereξ comes from the QNMs analysis in the right panel of Fig.3. We see that the red dashed line fits very well with the blue fitting line except in the large regime of m, such as m > 1.5. 2 This discrepancy originates from the differences between the two fitted values of b, i.e, b ≈ 0.1496 from Eq.(4.6) and b ≈ 0.1585 from Eq.(4.7). However, the error between these two b's is small, therefore, it is fair to say that our assumption of the relation between the coherence length and the graviton mass, viz. Eq.(4.6) is reliable. 2 In [50,61], the authors stated that the relation τ −1 s = γm 2 should hold as graviton mass is small, such that the momentum conservation is violated in a minor way. Therefore, if graviton mass is large enough, the above relation τ −1 s = γm 2 may not hold again. Consequently, the analytic relation in Fig.4 may not apply as m is large enough.

Conclusions and discussions
In summary, we have investigated the effects of momentum dissipation induced by the graviton mass on the critical dynamics in the background of a (3 + 1)-dimensional dRGT massive gravity by taking advantage of the AdS/CFT correspondence. The KZM scaling laws of the vortex number and the "freeze-out" time with respect to the quench rate were not affected by the graviton mass. This in turn indicated that the static and dynamic critical exponents were not changed by the existence of graviton mass. It further implied that the boundary field theory was always like a mean-field theory in spite of the graviton mass.
However, the outcome of vortex numbers depended on the graviton mass, since momentum dissipation would shrink the size of the symmetry-breaking domains. Inspired from Pippard's formula on the coherence length and the mean free path, we proposed an analytical relation between the coherence length and the graviton mass. This relation was verified from the QNMs analysis of the perturbations of the scalar field, which was a key result in our paper. We further investigated the relation between the number of vortices and the graviton mass both analytically and numerically, and the results were consistent with each other. Therefore, in conclusion, the momentum dissipation induced by graviton mass would increase the number of vortices, although it would not affect the KZM scaling laws.

A A Brief Review of dRGT massive Gravity
Let us consider the action for an (n + 2)-dimensional ghost-free dRGT massive gravity [34,57] where L is the radius of the AdS n+2 spacetime, m stands for the graviton mass parameter, f is a fixed symmetric tensor usually called the reference metric, c i are constants, and U i are symmetric polynomials of the eigenvalues of the (n+2)×(n+2) matric K µ ν ≡ √ g µα f αν : The square root in K means ( √ A) µ ν ( √ A) ν λ = A µ λ and [K] = K µ µ = g µα f αµ . The equations of motion turn out to be Rg µν − n(n + 1) One can choose a special form of reference metric as where h ij is the metric on a two-dimensional constant curvature space. A black hole solution of (3+1)-dimensional dRGT massive gravity with this reference metric is [69] ds 2 = −r 2 f (r)dt 2 + dr 2 r 2 f (r) + r 2 h ij dx i dx j , in which r + is the horizon of the black hole. In order to study the dynamics in this background, it is convenient to work under the Vaidya-AdS solution by making the transformation dυ = dt + dr/ r 2 f (r) ds 2 = −r 2 f (r)dυ 2 + 2dυdr + r 2 h ij dx i dx j . (A.10) By setting υ → t, r → 1/z, h ij → δ ij , k → 0, q → 0, (A.11) one can obtain the line element (3.1) in the main text.